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Preface 


Nuclear thermal hydraulics is the application of thermofluid mechanics within the nuclear indus- 
try. Thermal hydraulic analysis is an important tool in addressing the global challenge to reduce 
the cost of advanced nuclear technologies. An improved predictive capability and understanding 
supports the development, optimisation and safety substantiation of nuclear power plants. 


This document is part of Nuclear Heat Transfer and Passive Cooling: Technical Volumes and Case 
Studies, a set of six technical volumes and four case studies providing information and guidance 
on aspects of nuclear thermal hydraulic analysis. This document set has been delivered by Frazer- 
Nash Consultancy, with support from a number of academic and industrial partners, as part of 
the UK Government Nuclear Innovation Programme: Advanced Reactor Design, funded by the 
Department for Business, Energy and Industrial Strategy (BEIS). 


Each technical volume outlines the technical challenges, latest analysis methods and future direc- 
tion for a specific area of nuclear thermal hydraulics. The case studies illustrate the use of a subset 
of these methods in representative nuclear industry examples. The document set is designed for 
technical users with some prior knowledge of thermofluid mechanics, who wish to know more about 
nuclear thermal hydraulics. 


The work promotes a consistent methodology for thermal hydraulic analysis of single-phase heat 
transfer and passive cooling, to inform the link between academic research and end-user needs, 
and to provide a high-quality, peer-reviewed document set suitable for use across the nuclear 
industry. 


The document set is not intended to be exhaustive or provide a set of standard engineering ‘guide- 
lines’ and it is strongly recommended that nuclear thermal hydraulic analyses are undertaken by 
Suitably Qualified and Experienced Personnel. 


The first edition of this document set has been authored by Frazer-Nash Consultancy, with the 
support of the individuals and organisations noted in each. Please acknowledge these documents 
in any work where they are used: 


Frazer-Nash Consultancy (2021) Nuclear Heat Transfer and Passive Cooling, 
Study C: Reactor Scale CFD for Decay Heat Removal in a Lead-cooled Fast Reactor. 
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Introduction 


This case study demonstrates the use of Computational Fluid Dynamics (CFD) analysis to sup- 
port a continuing development process for the design of anew Advanced Modular Reactor (AMR). 
The overarching aim of this case study is to demonstrate the thinking and process for the suc- 
cessful planning, execution and utilisation of a thermal hydraulic analysis under typical industrial 
conditions. 


As with all the case studies in this series, this analysis provides a ‘worked example’ of a specific 
modelling task to illustrate the modelling approaches described in the technical volumes. Therefore, 
it is recommended that this case study is read in conjunction with the following technical volumes: 


Volume 2: Convection, Radiation and Conjugate Heat Transfer 
Volume 3: Natural Convection and Passive Cooling 

Volume 4: Confidence and Uncertainty 

Volume 5: Liquid Metal Thermal Hydraulics 


Many publicly available examples of thermal hydraulic analyses are performed as part of bench- 
mark studies (e.g. Study A: Liquid Metal CFD Modelling of the TALL-3D Test Facility). In these 
cases the geometry and boundary conditions are usually simple and well defined and validation 
data (often in the form of both high-fidelity analysis results and experimental measurements) are 
available. However, for many industrial thermal hydraulic modelling engineers, this does not gen- 
erally represent the situations in which they often need to work. 


Industrial simulations are built using the best information and methods available at the time, whilst 
recognising that there are unknowns and uncertainties that remain. There are often a number of 
non-technical constraints on the work, such as limitations on the time, budget and other resources 
available and there is a need to plan the analysis to align with these as well as meeting the technical 
objectives. 


This case study demonstrates the use of thermal hydraulic analysis within a ‘typical’ industrial 
scenario, including consideration of a range of technical and non-technical constraints. The specific 
objectives of this case study are to: 


¢ Demonstrate the planning and execution of thermal hydraulic analysis to add value to early 
stage reactor design. 


* Demonstrate the use of CFD to supplement information from system codes in a pool type 
Lead-cooled Fast Reactor (LFR) design. 


* Demonstrate how simplification can enable cost effective CFD modelling of the whole primary 
circuit of an AMR. 
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* Demonstrate how modelling assumptions and a lack of validation data affect the way that 
analysis predictions can be used. 


While the methods demonstrated by this case study use an LFR as an example, the specific details 
of the model and its results are not intended as important outputs and the intention is to assist the 
reader in developing their own models. Many advanced high temperature designs include pools or 
plena and decay heat removal systems that could be analysed in a similar manner and the stages 
needed to develop the analysis are similar to those required for any complex CFD model. 


To provide realism and improve the utility of this case study it has been developed in collaboration 
with Westinghouse Electric Company, who are currently developing a Nuclear Power Plant (NPP) 
based on LFR technology. As such, the layout of the reactor is broadly in line with the Westinghouse 
LFR design (as depicted by Ferroni et al., 2019). However, to enable open publication, various 
aspects of the geometry, component performance and analysis scenario are fictitious. Therefore, 
no quantitative data presented in this study (either input data or analysis results) should be taken 
as illustrating the actual design or performance of any previous, current or future Westinghouse 
LFR. 


The remainder of this document is structured as follows: 


Section 2 describes the LFR that is used in this case study, and the operating condition that 
is analysed using the CFD model. 


Section 3 describes the planning procedure for the analysis. This includes an identification 
of the results of importance, an assessment of the available input data and a selection of 
appropriate modelling tools. 


Section 4 describes the details of the analysis procedure. The modelling techniques applied 
to specific components in the LFR are described first, followed by the numerical methods 
that are applied in the CFD calculation. 


Section 5 describes the results that are extracted from the computational model, along with 
interrogation of those results in the context of the modelling assumptions and limitations. 


* Section 6 summarises the conclusions from the analysis in the context of a continuing reactor 
development process. Potential improvements to the model itself are highlighted, along with 
suggestions for additional analyses that could be carried out to provide further value. 
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2.1 


Defining the Problem 


The problem definition step in performing a thermal hydraulic analysis has a determining impact 
on its success and efficiency. A good understanding of why the analysis is being performed and 
what the results will be used for is vital to enable ‘fit-for-purpose’ decisions to be made throughout 
the model development. 


The scenario considered for this case study is an industrial setting at a specific point in time of 
the design of a new LFR. The target NPP output, basic layout and key operational performance 
parameters have been defined. The reactor core design is reasonably well advanced and initial 
estimates of the required performance of other reactor key components (e.g. the primary circuit 
pumps and heat exchangers) have been specified. 


The team have recently turned their attention to reactor performance under possible fault conditions 
and started to design a Passive Decay Heat Removal System (PHRS) to cool the reactor when the 
primary heat exchangers and pumps are not operating. Adequately demonstrating the capability 
of the system will be very important for the safety substantiation of the reactor and needs early 
consideration. 


The design of the reactor and PHRS are described in Section 2.1 and the specific industrial sce- 
nario and the need for the analysis are elaborated on in Section 2.2. 


Reactor Description 


The LFR considered by this case study is a pool type reactor, where the core and all primary 
components are contained within a Reactor Vessel (RV). A sketch of the LFR and the primary 
circuit flow paths is shown in Figure 2.1. 


Under normal operating conditions, the coolant passes up through the active reactor core and is 
heated by the fuel assemblies. The coolant exits the top of the core into the hot pool, where it 
moves upwards and radially outwards through the primary heat exchangers into the upper cold 
pool. Having passed through the heat exchangers, the coolant is then pumped into the lower cold 
pool and, from there, re-enters the bottom of the reactor core. 


The RV is surrounded by a Guard Vessel (GV). The GV is intended to retain the primary coolant in 
the unlikely event of a RV failure or a leak. Under these circumstances the lead would flow into the 
gap between the vessels and the width of the RV-GV gap is sized to ensure that the resulting lead 
level drop would not impair core cooling. 
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Water Pool 


ul 


Figure 2.1: A Sketch of the Pool Type LFR. 


The LFR considered in this case study is assumed to have the following dimensions. 


¢ The core barrel has an outer diameter of 4.0m. 

¢ The RV has an outer diameter of 8.0 m. 

¢ The GV has an outer diameter of 8.6 m. 

* The GV has an overall height of 10.1 m. 

* The active region of the core has a height of 1.25 m. 


In common with many AMR designs, this LFR incorporates a passive cooling system that can 
remove decay heat from the core without power or human intervention in the event of a fault. For 
example, in the event of a station blackout it is assumed that there is a loss of primary circuit forced 
circulation and the primary heat exchanger is no longer able to remove heat from the primary 
circuit. Natural circulation of the primary coolant (lead in this case) transports the decay heat out 
of the reactor core and, as the primary coolant flows downwards through the upper and then lower 
cold pool, radially outwards towards the RV walls. 


In this PHRS design, the decay heat is then transported across an argon-filled gap between the 
RV and the GV by convection and radiation. Finally, the decay heat is transferred from the GV toa 
volume of water on the outside of the GV by natural convection and boiling. If the fault persists for 
an extended period then the surrounding water will boil away. Natural convection of air then cools 
the GV outer surface and heat is transferred via multiple air chimneys to the environment. 


The design of the PHRS is important for both safety and the economic performance of the reactor: 


* It must be able to remove the decay heat during a bounding fault (defined by the safety case), 
without the RV overheating and compromising its structural integrity. 


* The PHRS design results in a continuous parasitic energy loss during normal operation that 
should be minimised. 
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2.2 


Defining the Problem 


These are competing requirements and so a satisfactory compromise must be reached if this 
design of PHRS is to be taken forward. 


Scenario Description 


Because the AMR under consideration is at a conceptual stage, some design decisions have not 
been made and several details are not refined. However, this does not mean that performing analy- 
sis at this point is without value. In fact, most modern reactors are first developed using knowledge 
from older reactor designs and modelling until the design reaches a level of maturity to justify the 
expense of specific experimental testing. Early analysis helps to improve understanding, gives con- 
fidence that the design will perform as intended, identifies areas of risk requiring further attention 
and feeds into the design of suitable test facilities. 


In the scenario postulated for this case study, system code modelling of the design under various 
fault scenarios has already been undertaken, which suggests that the proposed PHRS design will 
be successful. However, some aspects of the flow are in complex three-dimensional regions, such 
as the pools, which cannot be fully captured by this modelling approach. 


Stratification in the three reactor pools is a particular concern for the effectiveness of the PHRS, as 
excessive stratification in specific areas of the primary circuit may cause a blockage in the natural 
circulation flow path, reducing the flow rate and hence the effectiveness of heat transfer between 
the reactor fuel and the PHRS heat sink. Additionally, stratification has the potential to subject 
structural components to large thermal gradients, which maybe particularly problematic if they are 
fluctuating (i.e. thermal striping, which could cause damage by thermal fatigue). This phenomenon 
is not predicted with a high level of confidence by the existing system code model of this reactor. 
In particular, it is perceived that there is a higher risk of stratification under circumstances where 
the pool of water outside the GV has started to boil down, so the side wall of the GV is exposed to 
a high heat transfer rate below the water line and a much lower one above the water line. 
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Figure 2.2: Representative variation in core power and water level in the pool surround- 
ing the GV over the duration of the postulated fault scenario (station blackout). 


In order to address this concern, the modelling needs to consider a postulated fault that is suffi- 
ciently severe to result in boiling of the PHRS water pool. The chosen fault scenario is a station 
blackout (power failure). Figure 2.2 shows the variation in reactor thermal power for this fault, along 
with system code predictions of the water level after the initial loss of power. The core power drops 
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rapidly to around 3% of normal operation within a few minutes and less than 1% within the first few 
hours after the initial loss of power. The core power then reduces more slowly, reaching around 
0.5% of normal operation 6 days after the initial loss of power. 


At first the water level reduces more slowly than the core power because the water pool has to be 
heated to its boiling point before it can begin to boil. Once the water starts to boil, the water level 
reduces and reaches 50% of the height of the reactor (z = 2.4m) after around 1 day. The water 
level then continues to reduce and reaches the bottom of the reactor (when all the water has boiled 
away) after around 4 days. 


The requirement for further analysis is to supplement the existing system code results by providing 
further information on the possible impact of lead pool stratification (and other 3D flow phenomena) 
on the passive heat removal. At this stage in the project, the high level objective is sufficiently broad 
to enable exploratory analysis. However, a number of specific aims are refined around this high 
level objective. These specific aims are described in Section 3.1. The results of the analysis will be 
used to inform the design of the PHRS, to steer the direction of test activities relating to the PHRS 
and potentially to enable allowance for stratification to be incorporated into a future version of the 
system code model. 


It should be noted that the system code modelling itself is not the focus of this case study, therefore, 
details of it are not included’. Instances where input data is taken from a system code model, or 
where the results of different modelling methods are compared, have been described with the 
intention that a user may perform equivalent activities with the codes and results available to them. 
Information regarding system codes suitable for modelling the primary circuit of a Liquid Metal- 
cooled Fast Reactor (LMFR) are given in Volume 5, Section 3.2.2. To model boiling water in the 
PHRS, a containment code, such as GOTHIC, may be suitable. 


Several lead-based test facilities have been and are being operated to support development of 
LFR technology. However only a limited number of LFR design-specific test facilities have been 
constructed and there is no bespoke validation data. Therefore, it needs to be recognised from 
the outset that achieving a quantitatively accurate result with a high level of confidence will not be 
possible as part of this case study (neither is there much benefit in attempting to do this as the 
reactor design will likely change many times before it is finalised). 


Despite there not being onerous requirements for accuracy and confidence in the results, it would 
be helpful to identify areas of uncertainty that are likely to affect the conclusions and recommend 
activities that could increase this confidence in the future. There will be an ongoing requirement 
for thermal hydraulic modelling in the project and so a secondary objective of this initial work is to 
investigate how any new models can be used to add value beyond their initial scope. 


Since this case study is intended to represent one step in a fast moving industrial design process, 
a number of non-technical constraints are introduced to increase the applicability to real situations. 
The cost (both computational and in engineer time) and timescales for model development and 
analysis are often of high importance in an industrial context and cannot be ignored. It is common 
practice (and perfectly reasonable) to make compromises for these reasons, as long as they are 
recognised and their potential impact considered. 


1 Similar modelling using the same system code is described by Hua et al. (2020). 
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This case study includes the following non-technical constraints: 


¢ The initial analysis results are needed within approximately 10 weeks of starting the work. 

¢ The majority of the work is to be undertaken by an engineer with a moderate level of experi- 
ence with thermal hydraulic modelling. 

« The engineer has access to a specific set of existing modelling software (described in context 
in Section 3). 

¢ The engineer has access to limited computing resources, representative of what would be 
found in an industry setting where thermal hydraulic modelling is regularly performed. 
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After gaining an understanding of the need for modelling, it is then necessary to consider how to 
meet this need. In almost all cases, there will be a number of options and it is important to consider 
the implications of each decision on the ability of the model to meet its objectives. Even when time 
is short, thorough planning can significantly reduce the risk of error and rework thereby resulting in 
a reduction in the overall cost and timescales of analysis work. 


The Phenomena Identification and Ranking Table (PIRT) process, as described in Volume 1, Sec- 
tion 4.2 is a systematic approach to identifying the phenomena that have a dominant effect on the 
performance of a particular component or plant system (a relevant example is described by Liao 
et al., 2021). A formal PIRT process can be a large undertaking, involving a panel of experts to 
rank phenomena according to their importance and the level of knowledge associated with them. 
The requirement for analysis work is often identified from the outputs of a PIRT. When planning 
the analysis itself it is often helpful to use a similar mindset and list of considerations to identify (for 
example) the key modelling inputs and outputs and which phenomena will be most important. 


The planning activities form a sequence, with each gathering information and enabling decisions 
to be made. The information and decisions then flow down into subsequent activities. 


¢ Initial high level planning and information gathering identifies key model inputs, outputs, phe- 
nomena and selects an appropriate tool. In addition to identifying data that will be needed 
for the analysis, this planning confirms the feasibility of the modelling task and enables the 
selection of an appropriately skilled engineer to perform the work. 


* The next stage is to define the modelling strategy. This includes defining more precisely the 
number and nature of the analyses that will be needed and how they will be undertaken. This 
stage enables the development of a project plan for the analysis work and helps to identify 
risks in its execution. For some tasks, this stage of the planning is needed to calculate a cost 
(both in terms of engineer time and computational resources) and timescale for the work 
to feed into overall programme planning activities. This case study demonstrates another 
common scenario where the modelling strategy is (to some extent) driven by pre-defined 
time and resource constraints. 


¢ Finally it is valuable to plan the approach to quality control at the outset, noting that perform- 
ing the work in a traceable and rigorous manner is as important to the overall quality of the 
analysis as checking the results of the calculations. 


The sections below illustrate the planning process for this case study. While many aspects of this 
planning process are specific to this case study, a similar process could be followed for a wide 
range of thermal hydraulic analyses. 
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3.1 


Planning the Analysis 


Results of Importance 


The broad aim of the analysis was described in Section 2.2. However, more precision is now 
needed regarding what this particular model is trying to predict (and therefore what it does not 
need to predict). A complex model that tries to predict everything at once is likely to be both overly 
expensive and unsuccessful. One consequence of this is that a model developed for one purpose 
may not be as successful if used for another because modelling decisions that have a minor impact 
on some results may have a much larger impact on others. It is a common mistake to re-use 
an existing model for a different purpose without considering this issue and all the underlying 
assumptions and simplifications should be reviewed before this is attempted. 


In this case, the scenario of interest is a station blackout fault and similar faults have already 
been analysed using a system code model. The objective of the new modelling is to investigate 
the impact of three-dimensional flow phenomena that are challenging to calculate (or may not be 
possible to calculate) with the system code model. The following specific results of importance 
have been identified: 


* Primary coolant temperature distribution within and between the pools. 
¢ Identification of areas of stagnation and/or stratification within the primary circuit. 
¢ Primary coolant natural circulation flow rate. 


¢ Heat transfer across the cavity between the RV and the GV. 


It is noted that neither the fuel cladding temperature, nor that of any core components, are included 
in this list. This is not because these aspects are not important to the reactor operation, rather that 
this particular model is not aiming to predict them. 


It is often helpful for the thermal hydraulic modelling engineer to consider both the current require- 
ment and possible future uses for the model they are developing, as it is likely that there are steps 
that can be taken now that will make further development easier. However, judgement is needed to 
consider the cost of increasing the model’s capability where future application is uncertain. In this 
case, it is considered likely that, if successful, the model will be further developed to both improve 
confidence in the predictions and increase its functionality. The following results of future interest 
have been highlighted as potential areas for future development: 


¢ Temperature and mass flow distribution of coolant within the reactor core. 

* Distribution of coolant exiting the primary heat exchanger, including possible flow impinge- 
ment on the RV wall. 

* Temperature of important structural components, such as the RV and GV. 

¢ Analysis of reactor normal operation. 

¢ Analysis of different fault scenarios. 


There is no intention at this stage to develop a model to achieve the results of future interest. 
However, the impact of decisions and assumptions on these is considered throughout the model 
development process with a view to minimising the effort required for these future activities. 
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3.3 


Planning the Analysis 


Thermal Hydraulic Phenomena of Importance 


Having defined the results of importance, the next step is to consider what thermal hydraulic phe- 
nomena need to be predicted to achieve these results. At the same time, any known challenges 
in predicting these phenomena should be highlighted both to aid modelling tool selection and to 
inform the future interpretation of the results. Consideration of reactor operation under the fault 
conditions of interest and the specific results the model is aiming to predict reveals that several 
thermal hydraulic phenomena are likely to be important, including: 


¢ Natural and mixed convective heat transfer at surfaces. 

¢ Buoyancy driven flow phenomena such as plumes and stratification (as described in Volume 
3). 

Jet flow and mixing. 

Pressure losses, including frictional losses on walls and ‘form’ losses along the flow path. 


Turbulent dissipation of momentum and heat. 


Conjugate heat transfer and conduction through solids. 


Radiative heat transfer. 


From a modelling perspective, some phenomena will have a more significant impact than others on 
the results of importance and the effect of some are easier to account for with simple approaches 
than others. It is therefore not necessarily a requirement that the model predicts all of these phe- 
nomena to a high level of accuracy. 


It is often not possible to determine which phenomena will be most important and a level of judge- 
ment (often based on previous experience) is valuable at this stage. However, it is noted that a 
number of these phenomena (i.e. stratification and buoyant plumes) are judged to be both impor- 
tant and require a three-dimensionally resolved approach. It is also noted that the variety of thermal 
hydraulic phenomena listed above is likely to present a challenge and suggests that a modelling 
tool incorporating a wide range of functionality would be preferred. 


Input Data 


The predictions of any model are highly dependent on the input data available. If important in- 
puts (i.e. inputs that have a significant influence on the results of importance) are not known to a 
high level of confidence then the analysis results will also be of low confidence, regardless of the 
modelling method used. 


The input information available and identification of any gaps in knowledge has a deciding influence 
on what can be achieved with further analysis. In some cases, it will become clear at this stage that 
further analysis will offer no real benefits until more (or better quality) input information is available. 


To determine which input data are important, the results of importance and the phenomena of 
importance should be considered. For the scenario in this case study, important input information 
will comprise three-dimensional geometry, solid components thermal properties, thermal hydraulic 
properties of the fluids and a quantitative definition of the fault of interest. 
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The following input information is available for this case study: 


¢ The overall layout of the reactor. 

* The geometry and material of the RV, GV, upper cold pool, lower cold pool, core barrel and 
supports. 

¢ The geometry and layout of the core and location of control rod drivelines. 

¢ The size and locations of the primary heat exchangers. 

¢ The ducting for and location of the primary coolant pumps. 

« Estimates of the pressure loss characteristics of the core, primary heat exchangers and 
pumps. 

¢ Material properties of the solid components and fluids (density, thermal conductivity, specific 
heat capacity and dynamic viscosity for the fluids). 

¢ Details of the fault scenario of interest including the transient variation in decay heat load. 

¢ Results from and input data to a verified system code analysis of the fault of interest. 


The level of confidence with which these data are known is considered to be high for the geo- 
metrical information (although subject to change as the design progresses) and moderate for the 
pressure losses, material properties and fault scenario. The last item in the list above is of partic- 
ular significance. As the purpose of this study is to supplement the existing analysis (rather than 
perform an independent check of its results), some level of consistency with this analysis is impor- 
tant. There is also the opportunity to reduce the time (and therefore the cost) of setting up any new 
analysis by reusing a portion of the input data from the system code model. 


Investigation has also revealed gaps in, or inability to access the information which, when consid- 
ered in the context of the model requirements, may compromise its success (the impact of this 
missing information is discussed further in Sections 4 and 5). The following significant information 
is missing: 


¢ The detailed geometry of heat exchanger. 

¢ The detailed geometry of the pump (impeller and/or guide vanes). 

¢ The detailed geometry of the PHRS heat sink. 

¢ The geometry and associated pressure losses for the core exit and above core structure. 

¢ The geometry and associated pressure losses for the core inlet and below core structures. 
¢ Experimental data to enable validation of the model predictions. 


In all cases this information is not known at the current time because the design is still in progress. 
However, as the design progresses more of this information will become available and so could be 
incorporated into the model in the future. 


The extent of the data available is judged to be sufficient to build a geometrically resolved model 
capable of adding value to the reactor design process. However, the areas of uncertainty due to the 
missing information should be clearly highlighted in the results and specific studies to investigate 
the sensitivity of the model predictions to the assumptions made should be considered. 
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Modelling Tool Selection 


The objective of this case study involves the investigation of 3D phenomena and system code anal- 
ysis has already been undertaken, so a more geometrically detailed, three-dimensional modelling 
approach is justified. CFD modelling is likely to be able to give reasonable insight into the degree 
of stratification and its potential impact on the decay heat removal. Additionally, as discussed in 
Section 3.3 of Volume 5 (Liquid Metal Thermal Hydraulics), CFD has a greater potential to add 
value to the design of a pool type AMR than more traditional Light Water Reactors (LWRs). Using 
CFD at this point to contribute to reactor design activities provides an opportunity to explore the 
expense and utility of using this method in the future. 


However, using CFD in this context is not without challenges. As a geometrically resolved tool, the 
lack of geometrical information in areas where the design is not complete will require assumptions 
to be made. Additionally, the results of importance encompass the entire primary circuit, including 
prediction of the overall circulation rate, implying that a large portion of the reactor will need to be 
modelled. In order to develop one or more CFD models for a substantial portion of the reactor that 
can be solved with a reasonable computational cost, the size of the computational mesh will be 
limited and aspects of the reactor must be simplified even where the geometry is known. Both of 
these considerations will impact the accuracy of the results. Additionally, even with these simplifi- 
cations, an analysis of the full transient fault is not likely to be possible within the constraints of this 
study (as discussed in Section 3.8). 


In this case study scenario it is assumed that the modelling engineer is familiar with the ANSYS 
CFX CFD modelling tool and has access to a license (this being the preferred tool of their organisa- 
tion). This code contains most of functionality that is required for a model of this type. For example, 
it can easily include Conjugate Heat Transfer (CHT) and has a number of options for including the 
influence of buoyancy. Use of a code with which the modelling engineer has previous experience 
also helps reduce the risk of errors and facilitates more rapid model development. However, as AN- 
SYS CFX is a general CFD code, it does not include many special provisions for the heat transfer 
modelling of low Prandtl number fluids and this will need to be taken into account. 


It is not unusual in industry for an engineer to be limited to the range of software available within 
their organisation but, if timescales allow, consideration of more than one code is recommended as 
different CFD codes have different modelling options that might offer benefits to some applications. 
In this case, the modelling engineer does have access to a number of different CAD packages and 
a range of mesh generation software (discussed further in Section 4.2). 


Previous Work 


To aid in the development of the CFD model and to provide information as to what might reasonably 
be expected of such a model, it is helpful to investigate the availability of similar work. Volume 5, 
Section 3.4.2.2 lists a number of publications describing CFD analysis of reactor pools that can 
be used to inform the current modelling work. Comparisons with similar studies should not be 
considered as a form of validation, but are useful to build confidence nonetheless. 


A number of references have been identified as containing studies which are similar to that currently 
envisaged, and the following points of relevance have been noted: 
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* CFD studies of LMFRs have been undertaken in the past and resulted in a reasonable level 
of agreement with experimental results. 


* Findings from similar analyses suggest that CFD models are less accurate in the early stages 
of the fault, immediately after the reactor is shut down, than later in the fault when natural 
circulation is established (e.g. Zwijsen et al., 2019). 


* Similar studies generally comprise a single model encompassing the entire primary circuit 
and, due to the size of the modelling domain, tend to use a Reynolds-Averaged Navier- 
Stokes (RANS) turbulence modelling approach with wall functions for the boundary layers 
(e.g. IAEA, 2014). 


¢ A standard RANS turbulence model with a modified (increased) turbulent Prandtl number to 
reduce the turbulent transport of heat is favoured by many studies (e.g. Roelofs, 2019). 


¢ The size of the computational mesh used in similar studies generally results in areas of the 
reactor with a high level of complexity, such as the core, being represented by a porous 
media approach (e.g. Stempniewicz et al., 2019). 


In addition to these studies, Study A (Liquid Metal CFD Modelling of the TALL-3D Test Facility) pro- 
vides valuable learning regarding the modelling of liquid metals relevant to the primary circuit of an 
LFR. Study B (Fuel Assembly CFD and UQ for a Molten Salt Reactor) considers the modelling of a 
specific reactor component (a fuel assembly) under natural circulation conditions. Although Study B 
is concerned with molten salt, it provides useful information regarding the simplified representation 
of complex components. 


Extent of Model 


All thermal hydraulic models represent a finite portion of a complete system and need boundary 
conditions applied to the exterior surfaces of the model to represent the rest of the system and 
the surrounding environment. Defining the extent of a model is always a compromise between 
providing the most complete picture of the system and minimising the model size and complexity. 


In this case, the main objectives include predicting the primary circuit natural circulation flow rate 
and the three-dimensional flow phenomena in the pools. It will therefore be necessary to model 
each of the pools explicitly and the location, momentum and temperature of the flow into and out of 
each pool will be important. In addition, natural circulation involves a two-way coupling between the 
velocity and temperature fields (as described in Volume 3). In this case the heat source is the core 
of the reactor and the main heat sink is outside of the GV. When taken together, the requirements 
imply that it will be necessary to include the entire primary circuit, the RV, GV and PHRS heat sink 
in the model. 


In this scenario the regions outside of the RV offer potential for simplification. There is no require- 
ment for this model to provide detailed predictions of the temperature or velocities within the PHRS 
heat sink itself. The only result of interest relating to this component is the overall removal of heat 
via the PHRS. In addition, it is unlikely that the cover gas region above the surface of the pools or 
movement of this surface has a significant influence on the results of importance within the level 
of accuracy required of this model. Considering these effects, the model extent has been defined 
based on the following reasoning: 
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¢ The rate of heat loss from the fuel (and therefore the bounding heat load to be removed 
by the PHRS sink) and the large surface area of the GV is expected to result in nucleate 
boiling within the PHRS heat sink below the water level. Good correlations for heat removal 
via nucleate boiling are available and these can be applied as a boundary condition on the 
outside of the GV. 


* Above the water level in the PHRS, the quantitative value of heat removal by naturally driven 
air and steam is more uncertain and will depend on the size and shape of the PHRS cavity 
and air chimneys. Details of this cavity are not yet available so an attempt to model the flow in 
this region is not likely to give an improvement in accuracy compared with using correlations. 


* The heat transfer across the RV-GV cavity is not easy to consider using simple boundary 
conditions and will have a significant influence on predictions of PHRS performance. It is 
therefore concluded that this cavity should be included in the model. 


In addition to transferring heat to the PHRS, the primary circuit will also lose heat from the 
top surfaces of the pools. It is likely that applying simple boundary conditions to account for 
this will lead to some inaccuracies in the model. However, details of the geometry of this area 
of the reactor are not available so any explicit representation of the cover gas region would 
increase the model complexity (likely to the point where it could not be created and analysed 
within the timescales of this case study) with little increase in accuracy. Therefore, the cover 
gas region above the lead pool has been excluded from the model. 


The extent of the model chosen for this analysis is shown in Figure 3.1 and includes the reactor 
primary circuit, the gas in the RV-GV cavity and the solid components (e.g. RV, GV, core barrel and 
dividing plate). The boundary conditions applied to represent the PHRS heat sink and the cover 
gas region of the system are described in Section 4.1. 
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Figure 3.1: A sketch to show the extent of the computational model of the LFR. 
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The modelling of the entire primary circuit of a reactor using CFD represents a considerable chal- 
lenge for industrial computational resources. However, the symmetry of the reactor geometry lends 
itself to a sector modelling approach. The main thermal hydraulic components in the LFR un- 
der consideration are circumferentially periodic, with six primary heat exchangers and six reactor 
coolant pumps. Limiting the model extent to a one-sixth sector (as shown in Figure 3.2) with ap- 
propriately applied boundary conditions will reduce the total size of the mesh and therefore enable 
the analyses to be run with more manageable computational resources. 


Primary Heat Exchanger 


Reactor Coolant Pump 


Hot Pool 


Control Rod Driveline 


Cold Pool 


Figure 3.2: Diagram to illustrate the 60° periodic section of the reactor. 


It is noted at this point that just because the geometry of a given cavity is symmetrical does not 
necessarily mean that the flow patterns will be. Volume 3, Section 3.2.3 highlights that this is a par- 
ticular risk when dealing with buoyancy driven flows. This is considered further in the specification 
of boundary conditions (Section 4.1) and in the interpretation of the results (Section 5). 


Modelling Strategy 


The four case studies described in this series of documents present a number of different modelling 
strategies. This study builds on Case Studies A and B and the strategies, modelling approaches 
and findings of those cases are taken into consideration here. 


Study A presents a preliminary analysis, comparing RANS predictions with published Direct Nu- 
merical Simulation (DNS) data and then more complex analysis comparing both RANS and Large 
Eddy Simulation (LES) approaches. The incorporation of high-fidelity turbulence modelling meth- 
ods into this approach provided both additional confidence in the results and insight into the appli- 
cation of RANS turbulence models to the modelling of liquid metals. 


Study B presents a set of sub-models aimed at developing an understanding of the relevant phe- 
nomena in specific components before attempting a larger model encompassing a sector of the 
complete fuel assembly. This approach has some clear advantages in terms of early identification 
of specific modelling challenges and enabling more informed user-decisions to the made about, for 
example, turoulence models or mesh resolution. 
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However, consideration of these approaches in the context of the current study reveals that it 
would be unfeasible to perform this analysis task in entirely the same way as either Study A or B. 
The large geometrical extent of the model and the high number of potentially important thermal 
hydraulic phenomena would require a large number of initial studies to be performed to investigate 
everything. Furthermore, it is not immediately obvious which phenomena will be of high importance 
or how they will interact with each other in this case, making it difficult to focus sub-models to best 
advantage. In addition, research revealed no published DNS data that was sufficiently applicable 
to this study to use as a benchmark in the open literature. Using an LES approach for the analysis 
of even a subset of the modelling domain depicted in Figure 3.1 would be too costly within the 
constraints of this case study and would not be an efficient use of resources at this early stage of 
the reactor design. 


The modelling strategy demonstrated by this case study is to build a single model of the entire 
domain using reduced fidelity methods to represent the more geometrically complex reactor com- 
ponents. This approach is aligned with the level of input data available and has the potential to 
provide analysis predictions that are of use to reactor system design (as well as future model de- 
velopment) within the required timescales. These initial analysis predictions (and the experience 
of building and running the model) can be used to recommend focussed future studies and/or the 
development of more detailed sub-models to facilitate future modelling improvements. Despite the 
approach being described as a ‘single model’, development of the model prior to the final analyses 
is still likely to involve a progressive process of building and testing different regions and/or different 
aspects of the flow physics. 


The following approach to geometrical simplification is planned: 


* The geometry of the RV, GV, core barrel and dividing plate (between the upper and lower 
cold pools) will be modelled explicitly. The model will include the solid material within these 
components thermally coupled to the adjacent fluid regions to enable CHT. 


¢ The location and approximate geometrical extent of the core will be accurately included. 
However the details of the core geometry (e.g. the fuel assemblies) will be represented using 
a porous media approach (as successfully demonstrated by Study B and a number of the 
publications noted in Section 3.5). 


The total heat source provided by the core will be accurately captured by the model, with an 
appropriate spatial distribution within the core. 


¢ The control rod driveline channels passing through the hot pool will be explicitly represented 
where possible, noting that some simplification will be necessary close to the core outlet, 
where the detailed geometry is not yet defined. 


¢ The primary heat exchangers will occupy a physically correct total space envelope, with the 
inlet and outlet planes in the correct location and the inlet and outlet flow in the correct 
direction. Details of the internal structure and flow channels will be represented using porous 
media. 


¢ The primary circuit pump ducting will be accurately located, but the effect of the pump im- 
peller and any guide vanes will need to be represented using porous media. 


The following approach to capturing the relevant thermal hydraulic physics is planned: 
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¢ Due to the importance of both heat transfer and buoyancy, the variation in the fluid thermo- 
physical properties with temperature needs to be accounted for. The conductivity (and, for an 
unsteady analysis, the density and specific heat capacity) of the solid steel components will 
also be important for the prediction of PHRS performance. Therefore, variation in the thermal 
properties of all of the materials in the model with temperature will be included. 


The likely importance of buoyancy influenced phenomena require that buoyancy forces are 
captured as accurately as possible by the model. The system code modelling results suggest 
that the average lead temperatures within the pools do not show a large variation. However, 
local variation within each pool could be much higher, as could the variation in the argon 
temperature within the RV-GV gap. As recommended by Volume 3 (Section 2.3.2) a full 
‘temperature dependent density’ approach to modelling buoyancy will be taken. 


¢ The geometrical extent of the model and the non-technical constraints for this work necessi- 
tate the use of the RANS (or Unsteady Reynolds-Averaged Navier-Stokes (URANS)) turbu- 
lence modelling approach and (likely) the use of ‘wall functions’ to represent the boundary 
layers. The work carried out examining the performance of different RANS turbulence mod- 
els within Case Studies A and B and the recommendations within Volume 5, Section 3.3.2 
will be taken into account in the specific modelling choices. 


Although the final outcome is intended to be a single model encompassing the entire domain and 
all the relevant physical phenomena, the complexity of the model suggests that some preliminary 
analysis should be carried out as part of the process of building and self-checking the model. The 
following initial analyses are planned as a minimum: 


¢ A model of the RV-GV cavity with simple (but representative) boundary conditions to check 
that the radiation modelling is working correctly and to investigate the mesh resolution within 
the cavity. 


¢ Analysis of the primary circuit lead under forced flow conditions (with buoyancy disabled) to 
check the correct implementation of the regions of porous media. 


The results of these initial analyses have not been reported in detail (see Study B for a more 
comprehensive example of the use of preliminary models) but the findings of these analyses are 
referred to in the detailed model development described in Section 4. 


Fault Analysis Strategy 


This case study considers a fault that evolves in time, both in terms of a reducing core decay 
heat load and variation in the water level within the PHRS heat sink. Attempting to transiently 
analyse the complete fault using CFD, however, would require the use of a coarse mesh and 
significant computing resources (beyond that available to most industrial engineers). Even if this 
were possible, it would not be appropriate for a first attempt at CFD modelling of this reactor. 
On the basis that a system code analysis that does model the transient evolution of similar faults 
already exists, and that the aim of this work is to supplement those existing predictions, there is 
the potential to consider an alternative strategy. 


For this initial case, it was decided to analyse a ‘snapshot’ in time at a specific set of conditions 
during the transient. That is, the model would be analysed assuming that the heat load and PHRS 
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water level are not changing with time. The ability of this approach to capture a good representation 
of the transient flow patterns depends on the rate of change of the boundary conditions compared 
with the timescales of the flow phenomena being analysed. If, for example, the flow patterns and 
temperature distribution can establish themselves in a time-frame within which the external con- 
ditions have changed by only a small percentage, then this approach can be virtually identical to 
carrying out a true transient analysis. In a case where the external conditions are changing (or have 
recently changed) much more rapidly than the response time of the system being modelled, then 
this approach is unsuitable. Simple calculations estimating, for example, the residence time of the 
fluid in important areas or the thermal mass of fluid volumes or solid components and comparing 
these with the timescales of variation in external conditions can help determine the suitability of 
the approach. 
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Figure 3.3: Variation in the core power, water level and a system code prediction of the 
primary circuit lead temperature over the duration of the postulated fault scenario (station 
blackout) and the time snapshot analysed in this case study. 


In the context of this case study, the overall response of the ‘system’ (i.e. the reactor primary circuit) 
to the fault of interest can be estimated most easily from system code analysis results. Figure 3.3 
shows that the temperature of the lead in the reactor pools is predicted to rise rapidly over the 
first few hours of the transient due to the loss of cooling provided by the primary heat exchanger. 
However, after approximately 5 hours the initial heating is over and subsequent variation in lead 
temperature occurs more slowly, driven by a combination of the reducing heat from the core and 
the rate of heat extraction via the PHRS. The ‘snapshot’ approach proposed for this study assumes 
that the model is in thermal equilibrium (i.e. the heat flux into the model balances the heat flux out 
of the model, in the manner of a steady-state). Therefore, this approach is a better approximation 
when the average lead temperature is approximately constant. 
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It has been noted from previous published work that there can be relatively rapid transient changes 
in flow patterns (including thermal stratification) shortly after the initiation of the fault. For example, 
IAEA (2014) describes a number of CFD analyses of natural convection in the upper plenum of 
the Monju Sodium-cooled Fast Reactor (SFR). The benchmark results show the movement of 
stratification through the plenum over a timescale of around 15 minutes and a ‘snapshot’ approach 
would not have been appropriate for capturing such flow features. Although the layout and cooling 
arrangements of the reactor considered here are very different from Monju, it is likely that similar 
rapid transient changes in flow patterns could occur within this reactor in the first few hours of the 
fault. Therefore the early stage of the fault when the lead temperatures are rising rapidly was not 
selected for the ‘snapshot’ analysed in this case study. 


As the purpose of the new modelling is to investigate the primary circuit flow patterns when the 
water in the PHRS heat sink is partially boiled down, it is important to note that lead temperatures 
and core power are changing relatively slowly for the majority of time when the water is boiling, 
making the ‘snapshot’ approach feasible in this case. However, the water level is changing more 
rapidly (over 7 m in 2 days) so it would be wise to avoid times in the boil-down that may cause a 
significant change in the primary circuit flow patterns (for example when the water level is close to 
the dividing plate). 


The point in the transient chosen for this initial analysis (indicated by the arrows in Figure 3.3) rep- 
resents conditions when the PHRS water level is approximately half way down the upper cold pool. 
If this analysis is successful, a second point in time when the PHRS water level is approximately 
half way down the lower cold pool would also be of interest. Using a ‘snapshot’ approach at these 
times is unlikely to produce exactly the same results as a full transient analysis of the entire fault 
scenario. However, the analysis should still be able to provide insight into the overall flow patterns 
and presence of any persistent areas of stratification. 


It is noted that, even if the modelling is carried out with fixed boundary conditions, an unsteady 
(transient) analysis method may be required because a number of the expected flow phenomena 
cannot be well captured by a steady-state approach. 


The following analyses are planned once the model construction and the initial testing is complete: 


¢ A steady-state analysis of the chosen ‘snapshot’ conditions. 


¢ An unsteady analysis (using the steady-state above as a starting point) to assess the effect 
of unsteady flow phenomena. 


* An initial sensitivity study investigating the influence of the turbulent Prandtl number on the 
model predictions. 


Quality Assurance and Validation 


As described in Volume 1 (Introduction to the Technical Volumes and Case Studies) and Volume 4 
(Confidence and Uncertainty), verification and validation of any analysis are a key part of increasing 
confidence in the results. 


Because this scenario represents the creation of a completely new model in a tool not previously 
applied to this reactor, the modelling engineer will need to build the geometry and mesh and make 
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all of the user choices in the model for the first time. This increases the risk of errors compared 
with a situation where an engineer can base their model on previous verified work, where inputs 
and outputs are available for comparison. 


The purpose of the analysis at this stage is not specifically for the assurance of nuclear safety. 
However, the outcome forms part of the evidence that an organisation may use to assess commer- 
cial risk and make investment decisions regarding the continued development of the PHRS design. 
It is therefore important that the work minimises the potential for error and any shortcomings are 
well understood. 


It needs to be recognised that some aspects of CFD best practice, such as achieving a mesh 
independent solution, are not likely to be possible within the constraints of the current scenario 
and this should be clearly acknowledged to any users of the model predictions. In this case the 
following verification activities were undertaken: 


* Careful self-checking of all aspects of the model build and setup. 

* Progressive testing of sub-sections of the model during its development. 

¢ High level independent review of the overall modelling strategy. 

¢ Detailed independent checking of the geometry, mesh quality, input data, solution parameters 
and convergence. 

¢ Interrogation of the results to check the validity of any assumptions made. 

* Interrogation of the results to check for the influence of artificial boundary conditions and 
inadequate mesh resolution. 

* Independent review of the interpretation of the results, conclusions and recommendations. 


For any complex CFD model it is recommended that the verification is carried out in a staged 
manner. For example, it is worth checking that the geometry is ‘as intended’ before building the 
mesh to ensure errors are identified at the earliest opportunity and rework is minimised. 


Although no specifically relevant validation data have been found for this case, in its absence a 
number of other confidence building activities are possible, several examples of which are demon- 
strated in this report: 


* Comparison of key results with those from the system code analysis (Section 5.4.1). 

* Comparison of specific aspects of the analysis results with expectations from simpler valida- 
tion/oenchmark cases. 

* Sensitivity studies to investigate the importance of uncertain inputs that may have a signifi- 
cant influence on the results of importance (Section 5.4.2). 

* Comparison of conclusions with those of similar published work. 
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With an analysis plan in place and the input information gathered, the model can start to be devel- 
oped. It should be noted that the decision making process does not stop at the end of the planning 
stage and further refinement of the model implementation is expected throughout the development 
process. 


The basic stages of developing any CFD model involves the creation of the geometry, the cre- 
ation of the computational mesh and then the specification of the relevant CFD model inputs (e.g. 
material property data and boundary conditions). However, this process is not entirely linear, as 
decisions about the way in which a particular component/region will be defined in the CFD model 
can affect decisions about the geometry or the mesh. In this case study, the model development 
and analysis has been undertaken in the following stages: 


* The method used to simplify each component is described in Section 4.1. 

* The creation of the computational mesh is described in Section 4.2. 

* The material properties and other aspects of setting up the CFD model are detailed in Sec- 
tions 4.3 and 4.4. 

¢ The procedure used to achieve a solution (and monitor its convergence) is described in 
Section 4.5. 


Component Simplification and Modelling 


In an initial model, the geometrical simplification process often removes any small components that 
are not expected to have a significant effect on the flow, and reduces the level of detail explicitly 
included for the larger components. Even for more detailed CFD analysis than is to be attempted 
here, where a greater proportion of the geometry is known, simplification is often required. 


In addition to geometry, inputs to the model are also needed to represent other thermal hydraulic 
characteristics of the reactor components. For example, the decay heat from the core needs to be 
included. 


Throughout this simplification process, it is important to be mindful of the purpose of the analysis 
and the potential impact of any simplifications or assumptions made. 
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Circumferential Periodicity 


The main thermal hydraulic components in the LFR under consideration are circumferentially pe- 
riodic, with six primary heat exchangers and six reactor coolant pumps. Due to the resources 
available to undertake this analysis work, a decision was taken to utilise this circumferential peri- 
odicity to reduce the size of the model and the number of cells in the mesh. A one-sixth segment 
of the reactor has been used in this case. 


As discussed in Volume 3 (Natural Convection and Passive Cooling), buoyancy influenced flow 
may exhibit flow patterns that are not strictly symmetric and/or periodic even in a symmetric or pe- 
riodic geometrical space. Therefore, the boundary conditions that are applied to the circumferential 
sides of the one-sixth segment model need to chosen carefully. For example, hot plumes that rise 
above the core in the hot pool may drift around the centreline of the reactor. Applying symmetry 
boundary conditions to the sides of the segment would constrain this movement and could result 
in an inaccurate calculation of the flow field. For this example, periodic boundary conditions are 
preferred, as they allow flow patterns to develop that are not symmetric between the edges of the 
segment model. However, the use of periodic boundaries still provides a degree of constraint to 
the flow field, particularly near the centreline of the reactor, where the edges of the segment model 
are close together. Therefore the use of a periodic CFD domain should still be viewed with caution 
and the possibility of incorrect predictions of flow features such as plumes and convection cells 
considered in the interpretation of the results. 


Whilst it is usual for the major components in a pool type reactor to be circumferentially periodic, 
the fuel assembly pattern may not necessarily follow the same periodicity. If this is the case, the 
use of periodic boundaries may cause difficulties with the modelling of core. If the fuel assembly 
pattern does not follow the periodicity of the other components, it may still be possible to create a 
periodic domain by placing the periodic boundaries along the inter-wrapper gap, rather than cutting 
directly through the fuel assemblies (creating a jagged boundary; there is no need for a periodic 
boundary to be a single straight line). These arrangements are compared in Figure 4.1. 


Periodic Boundary 


(a) Straight Periodic Boundary (b) Jagged Periodic Boundary 
Figure 4.1: Different approaches for creating periodic boundaries for the core fuel assemblies. 


Due to the hexagonal shape of the fuel assemblies and the simple porous media approach intended 
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to be used to model the assemblies themselves, it has been possible to retain rotational periodicity 
with a straight slice along the boundaries (Figure 4.1a). All of the cut fuel assemblies are exactly 
one half except for the centre assembly and therefore connect ‘cleanly’ (in an easy to understand 
manner) to the opposite periodic boundary. 


The Core 


The core of a reactor is often the part that warrants the most detailed modelling effort. However, 
the focus of this case study is the investigation of stratification in the reactor pools during PHRS 
operation. Therefore, a high level of resolution of the internal structure of the reactor core is not 
essential and it can modelled as a source of heat and pressure loss. 


For this case study, a porous media approach to modelling the fuel assemblies is both appropriate 
and pragmatic. Porous media modelling is a common method in CFD, available in a wide range 
of codes. Further information regarding its use in reactor CFD analysis is described in Volume 2 
(Convection, Radiation and Conjugate Heat Transfer). Two levels of resolution were considered for 
modelling the core: 


* Concentric rings of porous material. 
¢ Hexagonal prisms of porous material to represent each fuel assembly. 


(a) Concentric Rings (b) Hexagonal Prisms 


Figure 4.2: Different core modelling approaches. 


These different core modelling approaches are shown in Figure 4.2. Both methods allow the heat 
load and flow resistance provided by the core to vary axially and radially and would be sufficient 
to provide consistency with the system code model. However, it was decided to use the hexagonal 
prism representation, for the following reasons: 


¢ With a hexagonal prism representation, a small gap can be created between each fuel as- 
sembly to allow inter-wrapper flow to be captured by the model. Learning from previous work 
such as Gerschenfeld (2019) suggests that the inter-wrapper flow can be important during 
passive decay heat removal operation, as it may provide convective transport of heat from 
the core to the core barrel. 
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¢ The distribution of the applied core heat load can be more accurately represented and, in the 
future, the core model could be developed to include different heat loads for specific fuel or 
control assemblies. 


¢ As further information becomes available regarding the core inlet nozzles and above core 
structure, it is likely to be easier to incorporate these into a model with a more realistic overall 
core structure. 


However, this method of representing the core has the following disadvantages when compared 
with the concentric rings approach: 


¢ Increased time and cost to generate a high quality mesh. 
¢ An increase in the total number of cells in the mesh. 
¢ It may be harder and/or take longer to apply some boundary conditions and input data. 


To allow the hexagonal prism representation of the core to incorporate the fuel assemblies and 
the inter-wrapper flow, it is necessary to consider how to divide the core into different regions 
without increasing the complexity of the CFD model beyond what is practicable. Some examples 
considering different porous ‘zones’ within fuel assemblies can be found in Yu et al. (2015), Wang 
et al. (2020) and Study B. In these examples, the assembly wrapper walls and inter-wrapper gap 
are modelled explicitly (i.e. without the use of porous media). However, both of these examples are 
concerned with predicting temperatures within the core itself and so may be over complicated for 
the purpose of this case study. 


Interwrapper Sandwich Porous Media Model 


Fuel Assembly Porous Media Model 
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Figure 4.3: Porous media modelling approach for fuel assembly and inter-wrapper gap. 


The inter-wrapper gap and the assembly wrapper wall are long and thin and it would require a 
large number of high aspect ratio cells to mesh them in detail. Some CFD codes offer the option of 
shell conduction to simulate heat conduction through and along the wall using thin (virtual) cells, 
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which are not explicitly resolved in the mesh. However, this treatment is not available in ANSYS 
CFX. An alternative option is to simplify these components in a similar manner to the fuel pins and 
consider them as a porous media. In this case study, it was decided to create a core model with 
two zones of porous media. The first zone includes the fuel pins and coolant inside the assembly 
and the second zone encompasses the fuel assembly wrapper and the inter-wrapper flow. This 
arrangement is shown in Figure 4.3. 


To allow the model to predict the distribution of flow between the fuel assemblies and inter-wrapper 
gap, different resistance coefficients were derived and applied to each of the porous regions. 


Fuel Assembly Pressure Loss 


For this case study, the process of specifying an appropriate resistance to the flow within the fuel 
assemblies can be facilitated by using the inputs already calculated for the existing system code 
analysis. 


The standard form of the loss coefficient defines the pressure loss as: 


ap =k (el?) 


where AP is the total pressure loss, p is the fluid density, U is the velocity magnitude. 


To specify consistent loss coefficients in the CFD code, it is necessary to understand how loss 
coefficients are defined. As different system codes and CFD tools do not use loss coefficients in a 
consistent manner, it is easy to introduce errors in the process. Taking information from the user 
manuals for both codes, the calculation in this case is presented below. 


This system code model uses lumped regions to represent the inner, middle and outer fuel as- 
semblies in the core, with loss coefficients specified per unit cross-sectional area squared of the 
lumped regions (giving overall units of m~*). Therefore, the loss coefficients from the system code 
were first multiplied by this cross-sectional area (Squared) to recover a standard form of the loss 
coefficient (K)'. 

K=¢xA? 


In CFX, pressure losses are applied in porous media as a negative momentum source term (mo- 
mentum sink) with units of force per unit volume. For example, in the vertical (z) direction, the 
momentum sink term is given by: 


1 
Se => ( 2 OU; + Kva30U2|UI) 


perm 
where Kperm is the permeability, Kis; is the inertial loss coefficient and yz is the dynamic viscosity. 


It is important to note that the velocity being used by CFX in this case is the so-called ‘superficial 
velocity’ which is the average velocity across the whole cross-section of the porous channel. In 
contrast, the local velocity through the reduced area of the real channel is referred to as the ‘true’ 


Some system codes use K = ¢ x 2A? so it important to check which form is applicable. 
This is consistent with the conventional use of pressure loss coefficients which often use the upstream mean or mid-channel 
velocity. 
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or ‘physical’ velocity and is not used here. Care should be taken to ensure that there is clarity about 
which one is being chosen and that the calculations are consistent. 


The pressure drop based on a standard loss coefficient is equivalent to a force per unit volume of: 


S,=- (*) (5ev-u') 


where L, is the length of the porous region in the vertical direction. It follows that a consistent 
pressure drop can be generated in the CFD model by applying the following value for the inertial 


loss coefficients: 
K 


L, 
Table 4.1 summarises the key values that were used in the calculation of K,.;, for the fuel assem- 
blies in different regions of the core. 


K loss = 


Region ¢ (m~“) A (m*) Ki) L (m) Koss (1/m) 
Inner Fuel Assemblies 15.2 1.96 58.1 2.47 23.5 
Middle Fuel Assemblies 18.9 1.80 61.1 2.47 24.7 
Outer Fuel Assemblies 8.61 3.45 102.4 2.47 41.4 


Table 4.1: A summary of the values that were used in the calculation of Kioss for the fuel 
assemblies in different regions of the core. A represents the cross-sectional area of the 
lumped region in the system code analysis. 


In the absence of further information it has been assumed that the control assemblies have the 
same porous parameters as the fuel assemblies. This is not likely to be entirely true in reality, but 
is judged to be a reasonable assumption for this initial analysis. 


In situations where pressure loss information is not available from a system code model, Study 
B (Fuel Assembly CFD and UQ for a Molten Salt Reactor) provides examples of deriving loss 
coefficients for porous zones in CFD models from both correlations (also used for the inter-wrapper 
gap below) and by the development of separate, detailed CFD models. 


Inter-wrapper Sandwich 


The inter-wrapper flow is not explicitly modelled in the system code analysis. Therefore, loss coef- 
ficients for the inter-wrapper flow cannot be extracted from the system code analysis and applied 
directly in the CFD model (as they were for the fuel assemblies). Instead, appropriate loss coeffi- 
cients for the inter-wrapper sandwich have to be derived from first principles. 


In the vertical direction, the inter-wrapper flow can be modelled as planar flow between two parallel 
plates. As this is a common type of flow, empirical data can be used to determine the pressure 
drop in the streamwise direction. This pressure drop can then be converted into an appropriate 
loss coefficient to apply in the CFD model. 


Before empirical data can be used, the flow regime needs to be classified as laminar or turbulent, 
so that an appropriate correlation can be selected. The system code analysis was used to provide 
an estimate of the flow rate through the core, and this flow rate was used to estimate the Reynolds 
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number (Fe) for the inter-wrapper gap: 


_ pUuD» 
UL 


Re 


where p is the density of lead, yw is the dynamic viscosity of lead, U is the cross-sectional average 
vertical velocity through the inter-wrapper gap and D,, is the hydraulic diameter of the inter-wrapper 
gap. As this calculation of the Reynolds number in the inter-wrapper gap is only an estimate, the 
calculation was checked after the predictions of an initial CFD simulation of the complete model 
were available. The Reynolds number in the inter-wrapper gap was found to be greater than 4000, 
and therefore the flow in the inter-wrapper gap can be considered fully turbulent. 


The pressure drop in the positive ‘z’ direction along the inter-wrapper gap can be calculated using 
the Darcy-Weisbach equation: 


boyd 
AP = —f (=) 5 Uz, true|U true| 


Dp 
where Uz true is the vertical-velocity averaged over the inter-wrapper flow passage area (the ‘true’ 
velocity), |Utue| is the magnitude of this velocity, L is the vertical length of the inter-wrapper flow 
region and f is the Darcy friction factor. As the flow is fully turbulent, there are a variety of empirical 
correlations that can be used to evaluate the Darcy friction factor. In this case study, the Colebrook- 
White correlation was chosen: 


1 — 910 ( € - 251 ) 
ie S10\ 37D, | RevF 
where e€ is the surface roughness of the fuel assembly walls (assumed to be zero in this case for a 
smooth wall). 


In the CFD model, the cross-sectional area of the inter-wrapper porous zone is larger than the inter- 
wrapper flow passage area, because it includes the thickness of the fuel assembly wrapper walls. 
Therefore, the superficial velocity calculated by the CFD code (U-) is lower than the true velocity 
through the inter-wrapper flow passage area (U2 tue). The superficial velocity can be written as: 


U; = Pa U,, true 


where gz is the area porosity of the inter-wrapper porous zone. The area porosity is the ratio of 
the flow passage area to the cross-sectional area of the porous zone and is less than 1. With this 
equation, the pressure drop through the inter-wrapper porous zone can be written in terms of the 
superficial velocity, so that it can be converted to a source term for the CFD code: 


i \ 2 L 
ap=—1(5,) seus (zr) 


This pressure drop is equivalent to a force per unit volume of: 


ae! 
S, = — (——x ) <pU,|U 
(sa) are 


It follows that a consistent pressure drop can be generated in the CFD model by applying the 
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following value for the inertial loss coefficient: 


f 


K joss er a ok 
Dnb4 


where the Darcy friction factor f is evaluated using the Colebrook-White equation. Table 4.2 sum- 
marises the values that were used in the calculation of Kj9<5. 


U (m/s) Dp (mm) Re (-) f (-) Pa Koss (1/m) 
0.1 Ta 5430 0.037 037 34.7 


Table 4.2: A summary of the values that were used in the calculation of Kis; for the inter-wrapper sandwich. 


If the flow through the inter-wrapper gap was found to be laminar (Re < 2500), then the pressure 
drop in the positive ‘z’ direction along the inter-wrapper gap can still be calculated using the Darcy- 
Weisbach equation: 


However, the Darcy friction factor is now calculated using the following equation: 


96 
f= — 
Re 


which can be derived using control-volume analysis for laminar flow between two parallel plates 
(Incropera et a/., 2011). By substituting this equation into the Darcy-Weisbach equation, it follows 
that the pressure drop in the positive ‘z’ direction along the inter-wrapper gap (if the flow was found 


48uL 
AP = — (Sar) U, 
Dida 


This pressure drop is equivalent to a force per unit volume of: 


48 u ) 
S22 (-0 WU: 
(oe 


It follows that a consistent pressure drop can be generated in the CFD model by applying the 


to be laminar) is : 


following value for the permeability in CFX: 


Dida 
48 


K perm — 


The walls of the inter-wrapper porous region were set to ‘free-slip’ walls (zero velocity gradient 
normal to the wall), rather than ‘no-slip’ walls (zero velocity). A ‘no-slip’ wall generates frictional 
losses, which are already accounted for by the momentum sink terms applied by the porous media 
model. The use of ‘no-slip’ walls in this case would therefore ‘double account’ for this frictional loss. 
The presence of a ‘no-slip’ wall also generates shear and (therefore) turbulence in a manner that 
may be inconsistent with the length scales of the component represented by the porous media. 
Turbulence in porous media Is discussed further in Volume 2. In this case study, ‘free-slip’ walls are 
applied to all porous zones in the CFD model. 
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Heat Transfer in the Core 


This CFD model does not need to (and is not able to) accurately predict the flow and temperature 
distribution within the core*. However, the representation of heat generation and heat transfer within 
the core should be sufficient to generate a reasonable prediction of the flow and temperature 
distribution entering the hot pool and the heat transfer to the core barrel. 


For simplicity, a thermal equilibrium modelling approach has been adopted for the solid and coolant 
within the fuel assembly porous zones and the inter-wrapper porous zones. With a thermal equi- 
librium model, the CFD code calculates a density, thermal conductivity and specific heat capacity 
for the porous medium, based on the volume-weighted average (volume-porosity) of the solid and 
fluid material properties. With this approach, the thermal mass of the porous zone is consistent 
with the thermal mass of the fuel assembly or inter-wrapper region. However, the solid and fluid 
components of the porous zone are calculated to be at the same temperature (i.e. this model does 
not calculate the temperature differences between the fuel and the local coolant). For this partic- 
ular model, this approximation is not likely to affect the results of interest. However, if this model 
were aiming to predict fuel cladding temperature or if the temperature differences between the fuel 
cladding and the coolant were expected to be high (for a less thermally conductive coolant for 
example) then this simplification would not be appropriate. 


The thermal conductivity of the fuel assembly and inter-wrapper porous media were both set as 
isotropic (equal in the vertical and lateral directions). It is possible to prescribe anisotropic prop- 
erties in CFD porous media. However, in this case, the superficial flow velocities in the core are 
expected to be low (= 0.1 m/s) and the coolant has a high thermal conductivity. Therefore the heat 
transfer along the length of the fuel assemblies is predicted to be of the same order of magnitude 
as the heat transfer across the fuel assembly (from the fuel pins to the wall). 


For coolants that are transparent (or partially transparent) to thermal radiation and have a signifi- 
cantly different conductivity to the solid walls, anisotropic conductivities are likely to be required for 
the porous media, even for an initial model. An example of such a coolant and how the anisotropic 
conductivities can be derived can be found in Study B (Fuel Assembly CFD and UQ for a Molten 
Salt Reactor). 


Core Inlet 


It is expected that components located below or at the inlet to the core will provide some resistance 
to natural circulation during a loss of core cooling event. Depending on the configuration, they may 
also provide a degree of flow turning and redistribution as the coolant enters the bottom of the core. 
As the effect of any below core structures may be significant, it would be advantageous to include 
a representation of them in the CFD model. However, at the time of this study, even approximate 
details of the core inlet geometry are not available and therefore cannot be included. 


Rather than neglect these structures completely, it was decided to include a block of porous ma- 
terial at the very bottom of the reactor core, as a placeholder (below the lower core in Figure 3.1), 
such that an appropriate pressure loss can be included in future versions of the model (when the 


This is due to the level of simplification used in the development of the core model and the lack of detail available regarding 
the above and below core structures (which will be key in determining the flow distribution). 
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relevant information is available). When the design has progressed significantly, it may be worth 
replacing this porous region with a more detailed representation of the geometry to capture further 
flow detail. In the current model this region is specified to provide no additional flow resistance. 


A consequence of the missing geometry at the inlet to the core is that the model predictions of 
flow distribution entering the core may not be representative of the true flow distribution. Hence, 
the flow and temperature distribution within the core may not be representative of the true distri- 
bution. However, as the primary purpose of this CFD model is not to assess the reactor core, the 
current modelling approach is considered to be acceptable providing that this limitation is taken 
into account in the interpretation of the results. 


Primary Heat Exchanger 


The primary heat exchanger in this reactor design is novel and many of the details are not available 
at the time of this analysis. It comprises an ‘active’ portion which removes heat from the primary 
circuit under normal operation and an upper body, which houses the secondary system coolant 
pipes (see Figure 2.1). The primary heat exchanger in this reactor design is located between the 
hot and upper cold pools. The flow moves laterally though the heat exchanger, exiting in a radially 
outwards direction towards the RV wall. 


During a loss of core cooling event, the primary heat exchangers are no longer the dominant source 
of heat removal from the primary circuit. In this analysis, a pessimistic approach was adopted and 
they were assumed to remove no heat from the reactor (requiring all the heat to be removed via 
the PHRS). Despite the lack of active heat removal, the primary heat exchangers are still large 
physical features in the upper cold pool which need to be included in the model. The tortuous flow 
path through the primary heat exchangers will also provide a resistance to the flow between the 
hot and cold reactor pools. 


Although detailed geometry of the heat exchanger is not available, the general design intent for this 
component, including: the approximate size, shape and location, the open area available for coolant 
flow and an estimate of the intended pressure drop through the heat exchanger are available. Given 
the available information and the need to simplify complex components where possible, the active 
portion of the heat exchanger has been modelled using a porous medium. 


The coolant leaves the heat exchanger radially outwards towards the RV wall. In reality, the flow at 
the outlet of the heat exchangers will form jets as the flow exits the individual passages within the 
heat exchanger (see Figure 4.4). In a porous media model the flow exits at the superficial velocity 
and will not resolve these jets. However, in this natural circulation case, the jets are expected to be 
at low velocity, losing momentum and mixing out shortly after leaving the heat exchanger. Hence, 
stratification in the cold pools is unlikely to be significantly affected by the jet structure at the heat 
exchanger outlet unless areas of high temperature gradient intersect with the heat exchanger outlet 
itself. Whist a simple porous media approach is appropriate for an initial model it could be improved 
in the future. Possible methods for future investigation include creating a more detailed sub-model 
of the heat exchanger and its local surroundings, or explicitly representing a portion of the heat 
exchanger close to its outlet, thereby allowing the flow to exit with the correct, physical velocity and 
direction. 
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(a) Outlet Ports Not Included (b) Outlet Ports Included 


Figure 4.4: Differences between including and not including a resolved heat exchanger 
outlet geometry in the CFD model. 


The resistance coefficients in the porous medium were specified to provide consistency with the 
assumptions in the system code analysis. The calculation of the porous loss coefficients for inclu- 
sion into the ANSYS CFX model is consistent with the method used for the core (as described in 
Section 4.1.2.1). In the vertical and lateral directions, the inertial loss coefficients are specified to 
be 10 times the magnitude of the inertial loss coefficient in the radial direction. In reality, the design 
of the heat exchanger may completely prevent flow in these directions or be quite open. However, 
this is considered a reasonable loss coefficient to apply at this stage and the overall flow is con- 
strained by the walls of the heat exchanger. It should also be noted that specifying an extremely 
high loss coefficient in any direction can lead to stability problems in the analysis convergence. 


The upper body of the heat exchanger contains the secondary coolant piping, structural supports 
and voids. The primary circuit lead can not enter this part of the heat exchanger. However, the 
upper body will still conduct heat and it will have thermal mass that is important for an unsteady 
analysis. 


In this case study, the upper body has been modelled as a steel box with a representative wall 
thickness. The inside of this box is filled with a block of lower conductivity material with a nominal 
specific heat capacity to represent the components contained within the body. This arrangement is 
shown in Figure 4.5. 


Reactor Coolant Pumps 


The reactor coolant pumps spool down within the first few minutes of a loss of power event. While 
the pump impeller no longer drives coolant around the reactor, the stationary impeller, internal 
geometry and pump housing still provide a resistance to the buoyancy driven flow that develops in 
the reactor. This resistance reduces the natural circulation flow rate and needs to be accounted for 
in the CFD model. 


In this case study, the detailed design of the pump is not available (although the approximate 
size, shape and location of the housing is known). Hence, it is not possible to create a detailed 
computational mesh of the internal geometry of the pump, to calculate the flow resistance or swirl 
it provides. Furthermore, as the design of the reactor develops, it is likely that the detailed pump 
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Figure 4.5: The modelling of the region above the heat exchanger. 


design will change to achieve the desired operational circulation rate. 


The resistance provided by the stationary impeller and pump internal geometry was included in the 
CFD model as a region of porous material (a porous-zone) within the pump housing, as shown in 
Figure 4.6. The pump housing has been included in the CFD model as a thin cylindrical wall. The 
level of detail of the geometry of the pump housing is necessary to ensure that the fluid enters the 
lower cold pool at the correct height, as the momentum of injected fluid could affect the stratification 
in the lower cold pool. 


The porous loss coefficients for the pump were specified to match the pressure loss characteristic 
used in the system code model. The loss coefficient in the vertical direction was derived using 
the same method that was used for the heat exchanger, and the loss coefficients in the radial 
and azimuthal directions were specified to be 10 times higher than the vertical loss coefficient (to 
promote model stability). In reality, the stationary impeller and pump internals may introduce swirl 
in the flow leaving the pump (i.e. locally, the loss coefficients will be lowest in a direction that is 
aligned with the impeller blades or vanes). However, it is not possible to estimate this effect without 
further details of the pump design. 


Porous-Zone 


Figure 4.6: Location of the porous-zone to represent impeller and pump internals. 
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RV-GV Cavity 


The argon-filled gap between the RV and GV forms an enclosed space and heat is transported 
across this cavity by natural convection and thermal radiation. For the current reactor design, the 
RV and GV walls are smooth and do not contain any fins or extended surfaces to increase the rate 
of heat transfer. 


The modes of heat transfer inside enclosures are discussed in more detail in Volume 2, Sec- 
tion 2.1.3. In this case study, two possible methods for modelling the heat transfer across the 
RV-GV cavity in the CFD model were considered: 


* Creating a separate fluid zone for the argon in the cavity and generating a mesh that is suffi- 
cient to resolve the convective flow patterns, thereby modelling the convection and radiation 
explicitly. 

* Modelling the cavity as a solid with an effective thermal conductivity to match the expected 
convective and radiative heat transfer. The required effective thermal conductivity can be 
calculated using an appropriate empirical or analytical correlation for the convection and 
analytical radiative heat transfer methods. 


The second method has a lower computational cost than the first. However, the first method en- 
ables more accurate modelling of both the convective and radiative heat transfer. Given the impor- 
tance of the prediction of heat transfer across the RV-GV gap and the potential for the temperature 
distribution of the RV wall to affect the stratification in the upper and lower cold pools, the first 
method was used in this case study. 


The system code results showed that the heat transfer across the RV-GV cavity is expected to be 
dominated by radiation. Therefore, an appropriate radiative heat transfer model needs to be en- 
abled in the CFD model. Because argon does not participate significantly in radiation, the Discrete 
Transfer radiation model (with the Surface-to-Surface (S2S) transfer mode) was selected in CFX. 


For this investigation, an emissivity of 0.9 was selected for both the RV and GV walls. This emis- 
sivity is representative of a painted surface finish and was chosen for consistency with the system 
code analysis. A sensitivity analysis for the effect of changes in the emissivity on the temperature 
distribution in the lead pool could be carried out with this CFD model. However, such a sensitivity 
analysis is not reported here. 


Outer Surface of the Guard Vessel 


During normal operation, the outer surface of the GV is surrounded by a pool of water. This pool of 
water acts as a sink for the heat transported to the GV wall. During a loss of core cooling event the 
primary heat exchangers are not operating, so the majority of the heat is removed from the outer 
surface of the GV by this pool of water. 


As described in Section 3.7, the point in the fault that is being analysed in this model is one 
where the pool of water has partially boiled down and part of the GV is exposed. This condition 
is of particular interest as the removal of decay heat below the water surface will be significantly 
higher than that above the water level, which may promote stratification. For this analysis, the water 
surface is located approximately half way down the upper cold pool. 
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The heat transfer from the outer surface of the GV is defined using boundary conditions in the CFD 
model. An external Heat Transfer Coefficient (HTC) boundary condition was applied to the outer 
surface of the GV. With this boundary condition, the heat flux (q) removed from the outer surface 
of the GV is given by: 

q = hext(Tw — Too) 


where hex; is the HTC on the outer surface of the GV, 7,, is the local temperature of the GV wall 
and 7. is a user defined sink temperature. 


The HTC and the sink temperature are specified by the user and remain constant throughout the 
calculation, while the wall temperature 7,, and the heat flux q are calculated by the CFD code and 
therefore vary during the calculation. The HTC can be specified to vary spatially with height, such 
that the heat transfer is significantly higher below the water level than above the water level. 


Above the water level, heat is removed from the GV wall by radiation and natural convection. The 
HTC for the natural convection of steam was calculated using an empirical correlation for a vertical 
flat plate (Incropera et al., 2011). 


2 
0.387Ra/® 
Nu = | 0.825 + “___. 
[1 + (0.492/Pr)9/16] 


Nuk 


Nese conv = L 


where L is the exposed height of the GV above the water level, k is the thermal conductivity, Pr is 
the Prandtl number and Ra is the Rayleigh number. The following formulation was used for Ra: 


gBb(Tw = Teal 
VQ 


Ra= 


where g is the acceleration due to gravity, 6 is the coefficient of thermal expansion, v is the kine- 
matic viscosity and a is the thermal diffusivity. k, Pr, B, v, w were all evaluated at the film tempera- 
ture on the outside surface of the GV wall (Tym = 0.5(Tw + Too))- 


The HTC for radiation above the water level was estimated by using the temperature of the GV 
from the system code results (7,,): 


eo (Tt — TS) 
ey ~~ Tx 


Next vad = 


where « is the emissivity of the outer surface of the GV and a is the Stefan-Boltzmann constant. 


In some CFD codes, it is possible to specify boundary conditions for convection and radiation occur- 
ring simultaneously (and this would be the preferred approach). In ANSYS CFX, this functionality is 
not available. Therefore, the heat transfer coefficients for radiation and convection were combined 
together to enable the CFD model to remove heat due to convection and radiation simultaneously. 


hext = Nee coi + Ae rad 


For this case study, the convective contribution was estimated to be approximately 5 W/m?K and 
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the radiative contribution was estimated to be approximately 40 W/m°K. Therefore, both the con- 
vective and radiative contributions were deemed to be essential to include in the thermal boundary 
conditions. 


Although this approach is considered sensible for the model at this stage, it is an area where the 
model could be improved in the future: 


* There will be a layer of droplet laden vapour on the GV wall above the water surface. Vapor- 
isation of the droplets will provide additional cooling to the guard vessel wall, increasing the 
overall heat transfer coefficient. 


The droplets will participate in the radiative heat transfer above the water surface. 


¢ As the water boils away, ambient air is drawn into the cavity to replace the steam which 
passes out of the cooling chimneys. The reactor decay heat is increasingly removed from 
the GV walls by natural circulation of air and radiative heat transfer. Therefore the relative 
contribution of the convective and radiative heat transfer mechanisms will need to be re- 
assessed as the water boils away. 


Below the water level, heat is removed by boiling the water pool. The HTC for boiling is orders of 
magnitude higher than for natural convection and radiation above the water level. It is therefore 
expected that the majority of the decay heat is removed from the GV below the water level. 


To specify a HTC below the water line, empirical correlations for boiling can be used. To use an 
empirical correlation, the mode of boiling has to be identified first. In this case study, the core decay 
heat was divided by the surface area of the GV below the water line to give an approximate heat 
flux (per unit area). This approximation neglects other heat sinks but is sufficient to identify the 
likely mode of boiling. It was concluded that the regime is likely to be nucleate boiling. 


In this case study, the correlation developed by Rohsenow for nucleate pool boiling (Incropera 
et al., 2011) was selected for simplicity. In this correlation, the temperature difference (Tw — Tsar) 
required to sustain a heat flux of q is given by: 


q=p h &(p! _ pv) ae Cail Ty _ T sat) 
eve oO Cophee Pr, 


where the subscript / is for the liquid phase and the subscript v is the vapour/gaseous phase. 
Table 4.3 summarises the quantities and the notation used in the Rohsenow correlation. 


In the CFD model, the boiling heat flux is applied to the outer surface of the GV as an HTC: 


With a temperature difference (Tw — Tsar) of (Say) 20°C between the wall and the water, the HTC 
below the water level was calculated to be approximately 50, 000 W/m°K. As expected, this is sev- 
eral orders of magnitude higher than above the water level, suggesting that the majority of the 
decay heat in the CFD model will be removed below the water level. It is noted that other correla- 
tions may give higher or lower values of HTC. However, the precise value is likely to be unimportant 
whist it remains several orders of magnitude greater than the effective HTC across the RV-GV gap. 
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Quantity Notation Units 
Latent Heat of Vaporisation heg J/kg 
Acceleration due to Gravity g m/s? 
Density p kg/m? 
Surface Tension Coefficient o N/m 
Specific Heat Capacity Cp J/kg K 
Prandtl Number Pr - 
Surface Temperature Tw K 
Saturation Temperature Tsat K 
Empirical Coefficient Cag - 
Empirical Coefficient n : 


Table 4.3: Notation used in the Rohsenow correlation. The empirical coefficients C.;,¢ 
and n have been taken as 0.0132 and 1.0 respectively, for water boiling on a stainless 
steel surface (Incropera et al., 2011). 


For the initial model development, the boundary condition was specified using a single constant 
value for each of hex; and 7,3:. However, a more sophisticated approach using the correlation above 
to apply a temperature dependant HTC could be incorporated in the future for increase accuracy. 


Free Surface of the Lead Pools 


After the initial loss of power event, the level of the lead pools will change, due to a decrease in flow 
rate (i.e. a decreased pressure difference between the hot and upper cold pools) and an increase 
in the bulk temperature of the lead (thermal expansion). 


In this analysis, an isolated snapshot in time is assessed, representing 0.5 days after the initial loss 
of power event. At this point in time, the system code analysis predicts that the natural circulation 
flow is established and the lead temperatures are changing very slowly. The depth of the lead 
pools will therefore be approximately constant and the low flow velocities suggest that agitation 
of the free surface is likely to be minimal. Therefore, the free surfaces of the lead pools were 
modelled as flat, rigid surfaces, with zero gradient (free-slip) boundary conditions applied normal 
to the surface in the momentum equations. This boundary condition is equivalent to setting zero 
shear stress between the lead and the argon at the free surface. 


Due to the large surface area, heat loss from the surfaces of the pools may make a significant 
contribution to the overall heat balance for the reactor. In reality, heat will be lost from the free 
surface by both convection to the argon cover gas and radiation to the structures above the sur- 
face of the lead pools. However, in both cases the magnitude of the heat loss will depend on the 
geometry, temperatures and layout of the cover gas region, which was not available at the time of 
this analysis. 


The lead in the reactor heats up considerably during the course of the fault, so the heat transfer 
from the surface of the pools to the area above is likely to be dominated by radiation. A very simple 
approach to modelling the thermal radiation is to assume that the area of the surfaces above the 
pool is large (the top cover plate above the lead pools will contain penetrations and there is likely 
to be a significant number of components passing through the cover gas region) and therefore 
they behave as a black body. This is not likely to be completely true in reality, but is a necessary 
assumption without further information about the components in the reactor cover gas region. The 


36 of 105 


4.2 


Performing the Analysis 


view factor from the pool surface to the surfaces above will be 1 since the surface cannot ‘self view’. 
The heat flux (per unit area) then becomes: 


ext,rad = O€ (ar a T) 


where « is the emissivity of lead surface, o is the Stefan-Boltzmann constant and T.. is a user- 
specified bulk temperature representing the solid components above the surface of the lead. 


The surface emissivity of the lead is difficult to specify as it is likely to be affected by the presence 
of lead oxide floating on the surface of the pool. Additionally, 7,, is difficult to specify in absence 
of any information regarding the equipment above the surface of the pool. For this initial case, 7. 
was specified as 200°C and « was specified as 0.95. 


In some CFD codes, heat loss boundary conditions due to radiation can be specified directly us- 
ing the emissivity of the surface and external thermal radiation temperature (T,.). However, CFX 
does not have this functionality. Therefore, the heat loss boundary condition was specified as a 
temperature dependent heat transfer coefficient instead: 


Next tad = 


The approach used in the case has intentionally been kept very simple and should certainly be 
revisited when further information becomes available for the region above the lead free surface. 
At this stage the effect of this boundary condition on the results of importance should considered 
in the analysis of the model results and a sensitivity study may be required to investigate the 
importance of the assumptions made. 


Meshing 


Before constructing a computational mesh, an appropriate level of resolution has to be chosen. In 
general terms, the level of resolution is chosen based on a balance between: 


* Increasing the number of cells in regions of interest, which increases the resolution and 
accuracy of flow features. 


¢ Reducing the total number of cells in the mesh, which reduces the computational cost of the 
simulation. 


Since the mesh resolution required depends on the flow features, it is not possible to completely 
specify its adequacy in advance of the analysis. It is common practice to start with a mesh that is 
potentially coarser than needed and, upon analysing the results, identify areas where further mesh 
refinement would be beneficial. 


Some examples of where this has been done and the judgements made in the development of 
the mesh are provided throughout this section. However, it is often possible to estimate a reason- 
able starting point from knowledge of the modelling parameters (e.g. the approach to modelling 
turbulence) and expected flow patterns. More details of specific requirements for CFD meshing are 
described in Volume 1 (Introduction to the Technical Volumes and Case Studies) and Volume 3 
(Natural Convection and Passive Cooling). 
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This case study is intended to be representative of a typical CFD analysis that may be carried 
out during the early stages of reactor development, with limited resources available. Therefore, a 
high priority has been given to selecting a level of mesh resolution that would facilitate an efficient 
analysis within the required time frame, while maintaining a reasonable level of resolution in the 
main areas of interest. It is intended that a RANS turbulence modelling approach will be used and 
wall functions will be employed in the boundary layers. 


The mesh generating software ANSA was chosen to generate the mesh for this case study. ANSA 
was chosen for the following reasons: 


* The engineer has access to the software and is confident with its use. 


¢ ANSA can generate and connect a variety of cell shapes, structures and patterns and gives 
the user a high level of control over the mesh. 


* The chosen CFD solver (ANSYS CFX) cannot operate on meshes that contain polyhedral 
cells. ANSA is able to generate meshes that are suitable for CFX, which are dominated by 
large areas of hexahedra joined entirely by pyramids and tetrahedra. 


* ANSA is well suited for creating conformal interfaces between the solid and fluid regions, 
which can reduce interpolation errors between zones in the mesh. 


¢ ANSA is able to generate periodic boundaries with matching faces. Matching faces across 
periodic boundaries reduces interpolation errors and is required in some CFD codes which 
cannot process periodic faces that are not identical. 


¢ ANSA can copy and duplicate different mesh regions, maintaining links between them. This 
is particularly useful for the fuel assembly porous zones, which are geometrically identical. 


The resulting total mesh size in this case is around 13 million cells, which can be solved using 
computing resources that are available to most engineers in a reasonable time frame. For example, 
with 40 parallel processing cores, a steady-state calculation with 5,000 iterations was completed 
in 6 to 10 hours of real time, while a transient calculation of 20 seconds of simulated time was 
completed within a few days of real time. This speed of computation is practical for most engineers 
and does not necessitate access to larger dedicated computing resources with hundreds of cores*. 


It was noted early in the meshing activities that CFX uses the computational mesh a little differently 
from some CFD codes, by converting it into control volumes centred on the nodes of the original 
mesh before performing the calculation (further information is provided in the code documentation). 
However, since this conversion results in a similar (or slightly more refined) level of resolution, no 
specific changes to meshing approach are required to accommodate this feature of the code. 


The remainder of this section details the mesh structure and resolution in different regions of the 
reactor. Some brief calculations are also presented to demonstrate how an appropriate level of 
mesh resolution was estimated in each area. 


The time needed to converge a complex CFD analysis is more complicated than only considering the number of cores 
available and the mesh size. Other aspects of the hardware, such as memory, storage and interconnect speed will also 
make a difference, as will the specific choice of code and the temporal and spatial characteristics of the phenomena within 
the solution. It is recommended that a model is initially run to assess the performance and the way in which convergence 
progresses with iterations, before a larger number of analyses are undertaken. 
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RV-GV Cavity 


The argon-filled cavity between the RV and the GV is modelled explicitly as a fluid region. The 
guidelines for constructing CFD meshes for natural-convection driven cavity flows are less well 
defined than for other CFD applications. This is partly because it is difficult to clearly define (or 
estimate in advance) the momentum and thermal boundary layer thicknesses within the cavity. 
The flow is also likely to be three-dimensional, particularly if convection rolls develop. Therefore, a 
series of simple requirements were used to develop the initial mesh: 


¢ The wall adjacent cell should be placed within the viscous and conduction sublayers; a target 
y* of less than 5 was used. 


¢ Use a sufficient number (and distribution) of cells across the thickness of the cavity to resolve 
the flow structures. 


¢ Use a sufficient number (and distribution) of cells along the length of the cavity to resolve the 
flow structures and maintain a reasonable aspect ratio. 


¢ Use a growth ratio normal to walls such that sudden/large changes in cell volume between 
adjacent cells is avoided. 


¢ Use hexahedral cells where possible, as tetrahedral cells can lead to artificial numerical 
diffusion of the gradients in flow variables. 


Before constructing the mesh, the likely flow conditions in the cavity were determining by calcu- 
lating the Rayleigh number (Fa) for the flow in the cavity. The following formulation was used for 
Ra: 


_ gBTH- Te)L? 
Va 


Ra 


where g is the acceleration due to gravity, B is the coefficient of thermal expansion, Ty is the 
temperature of the outer surface of the RV, Tc is the temperature of the inner surface of the GV, L 
is the thickness of the cavity, v is the kinematic viscosity and a is the thermal diffusivity. B, vy and a 
were evaluated at the average of the RV and GV wall temperatures (0.5(Ty + Tc)). 


Ra was calculated to be around 2 x 10°, indicating that turbulent flow is likely to be present in the 
cavity. For some unconfined flows (for example, forced convection over a flat plate), an estimate of 
the required first cell height to achieve yt ~ 1 can be made using empirical correlations. However, 
such empirical correlations are less applicable for the enclosed buoyancy-driven flow that occurs 
in the RV-GV cavity. Therefore a different, geometrical approach was used to estimate the mesh 
resolution required in this case. 


Due to the narrow width of the RV-GV cavity and low expected flow velocities, it is likely that the 
boundary layers on the walls will occupy much of the volume of the cavity. For CFD applications, it is 
usually recommended that at least 10 to 15 cells are applied across the thickness of the boundary 
layer. Using this recommendation as a guide, 20 cells were specified across the half-width of the 
cavity. A geometric growth ratio of 1.2 (chosen based on the past experience) was then applied 
normal to the RV and GV walls, to reduce the height of the cells close to the wall. With these 
parameters, the height of the wall adjacent cell was calculated to be around 1mm. The resulting 
mesh in the RV-GV cavity is shown in Figure 4.7. 


A preliminary CFD analysis was carried out to assess the resolution of the flow in the cavity with 
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the chosen mesh parameters. The results of this preliminary analysis showed that the chosen 
level of mesh resolution was sufficient, with yt ~ 2 and that the flow structures within the cavity 
were sufficiently resolved. At this stage, if y* was found to be larger than 5, or flow structures 
within the cavity were not sufficiently resolved, then the mesh in the RV-GV cavity could have been 
refined, before running the final model. The flow in the RV-GV cavity predicted by the final model 
is examined later in Section 5.2.1. 


Reactor Vessel Guard Vessel 


CC 
40 cells 


Figure 4.7: Close-up view of the mesh in the RV-GV cavity. 


Lead Pools 


Mesh Structure Away from the Walls 


It is difficult to predict the flow structures in the lead pools that require increased resolution before 
carrying out the CFD analysis. Therefore, it is difficult to develop a bespoke mesh with targeted 
refinement in specific areas, and so a simple meshing approach was adopted instead. A relatively 
uniform distribution of cells was generated throughout the bulk of the pools, with increased resolu- 
tion close to the walls for the anticipated momentum and thermal boundary layers. Throughout the 
mesh, large changes in adjacent cell sizes were avoided. 


To generate a relatively uniform distribution of cells throughout the bulk of the lead pools, a ‘hexcore’ 
approach was used, with low aspect ratio 25mm hexahedral cells, as shown in Figure 4.8. A 
hexcore approach was chosen for the following reasons: 


* Hexahedral cells reduce the numerical dissipation of the thermal gradients (when compared 
to tetrahedral cells). Hence, they are more likely to capture stratification in the upper and 
lower cold pools. 


¢ The (approximately) cube shaped hexahedral cells have good mesh quality metrics. This 
improves the stability of the equations solved by the CFD solver. 


The ability to control the mesh generation to the extent needed to build this structure was one of 
the key drivers behind the choice of mesh generation software chosen in this case study (ANSA). 
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Inflation Layers 
Control Rod Hot Pool 


Fuel Assemblies 


Inflation Layers 


Figure 4.8: Close-up views of the mesh in the lower cold pool and the hot pool. 


An unavoidable consequence of the hexcore meshing approach is that a small ‘buffer’ region is 
required between the inflation layers (prismatic cells) that are grown from the wall to accommodate 
the momentum and thermal boundary layers, and the hexahedral cells away from the walls. The 
buffer region is filled with tetrahedral and pyramid cells, to allow the mesh structure to transition 
from inflation layers to the hexcore. While this buffer region does lead to a local reduction in cell 
quality, there were no noticeable artefacts in the flow and temperature fields in the lead pools in 
the vicinity of the buffer region. Therefore, the hexcore meshing approach was deemed to be the 
best option for this case study. 


Mesh Structure Close to the Walls 


The mesh structure close to the walls was generated with the aim of increasing the resolution of 
the momentum and thermal boundary layers on the walls of the lead pools. Inflation layers were 
grown using a geometric series, characterised by a First Cell Height (FCH), geometric growth ratio 
and number of layers. The FCH was specified first, with the aim of placing the first cell in a specific 
region of the momentum and thermal boundary layers. As discussed in Volume 5 (Liquid Metal 
Thermal Hydraulics), in a low Prandtl number fluid such as lead, the thermal boundary layer is 
considerably thicker than the momentum boundary layer. Therefore, the height of the first cell ina 
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low Prandtl number fluid corresponds to a different region in the momentum and thermal boundary 
layers. 


To resolve the boundary layer with the computational mesh, it is necessary to place the first 
cell within the viscous (and conduction) sublayer. In some analyses which use CFD (aerofoils 
at high angles of attack for example), it is necessary to do this and the FCH is specified to target 
yt ~ 1. However, when the modelling domain includes many solid surfaces, this approach greatly 
increases the number of cells in the mesh. The use of wall functions to model regions of the bound- 
ary layer close to the wall can significantly reduce the number of cells in the mesh, at the expense 
of a reduction in accuracy. 


In this case study, the FCH was specified to target y* in the range of 30 to 100 (ideally 30 to 60 
on the majority of the surface), which corresponds with the logarithmic region of the momentum 
boundary layer and the conduction sub-layer of the thermal boundary layer. This level of resolu- 
tion is a compromise between the overall time taken to complete the CFD analysis and accuracy 
required at this stage in the reactor development. Later in this case study, once the CFD analysis 
has been carried out, the flow and temperature profiles are examined in the vicinity of the wall to 
assess whether further refinement would be beneficial. 


It is difficult to estimate the required height of the first cell in the lead pools to achieve a target y™, 
since the flow in the lead pools is complex and it is not clear which empirical formulae would be 
appropriate. Therefore a more direct method was used. A preliminary CFD analysis was carried 
out with an FCH of 5mm (chosen using engineering judgement). The observed values of y* were 
in the range of approximately 150 to 350. The first cell height was then refined from 5mm to 1.5mm 
to achieve y* values in the target range. 


Having specified the height of the first inflation layer, a geometric growth ratio of 1.2 was selected 
for the remaining inflation layers. A value of 1.2 was chosen based on past experience of the 
analysis of internal flow systems. The number of inflation layers were chosen by trial-and-error, 
with the aim of minimising the transition in cell volume between the final inflation layer and the 
hexcore mesh. 12 inflation layers were found to be sufficient to minimise the cell volume transition 
in this analysis. 


Core 


The flow through the reactor core is one of the most complex regions of the reactor and would 
require a fine mesh to accurately capture the flow features. However, the purpose of the CFD model 
in this case study is not to resolve the flow features within the core. Therefore, the mesh only needs 
to be of sufficient resolution for the application of the heat load and to calculate the pressure drop 
across the core, so that a reasonable representation of the flow entering and leaving the reactor 
pools is produced. 


As discussed in Section 4.1.2.1, the core model was divided into regions of porous material to rep- 
resent the individual fuel assemblies and the inter-wrapper gap. This approach allows the different 
heat loads and loss coefficients to be applied to each fuel assembly individually, and allows some 
Due to the large thickness of the thermal boundary layer in low Pr fluids, it is likely that it will extend beyond the wall 


prismatic cells in some places. This will not necessarily cause any problems with the solution, but interrogation to check for 
any artefacts in the results caused by the change in mesh structure is recommended in areas of high importance. 
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representation of the inter-wrapper flow to develop in the CFD model. Figure 4.9 shows a close-up 
view of the mesh in the core to highlight the different regions of porous material. 


As the flow through the porous zones in the core (especially within the fuel assemblies) is expected 
to be predominantly in the vertical direction (along the axis of the fuel assemblies), a swept meshing 
approach was used to generate high quality prismatic cells, which are aligned with the vertical 
direction. Each fuel assembly was meshed individually in 2D first, and then swept in the vertical 
direction to generate the 3D mesh. 


Shield Assemblies 
(Yellow) 


Interwrapper 
(Green) 


Bypass 
Region 


Periodic 
Boundaries 


Fuel Assemblies 
(Orange) 


Interwrapper 
(Green) 


Slip Wall 


Figure 4.9: Close-up views of the fuel assembly and bypass region mesh structure. 


The inter-wrapper porous zone was generated using the same approach as the fuel assemblies. 
A 2D mesh was generated first and then swept in the vertical direction to generate a high quality 
prismatic mesh. As shown in Figure 4.9, the inter-wrapper porous zone was constructed to be 
fully conformal with the surrounding fuel assemblies and also the lead pools above and below the 
core (Figure 4.8). This conformal approach was adopted to minimise interpolation errors in the 
velocity field between the leads pools and the inter-wrapper mesh, as the inter-wrapper cells are 
considerably thinner than the hexcore cells in the lead pools. 


The inter-wrapper porous zone was meshed with three fluid cells across the thickness between 
adjacent fuel assemblies (CFX converts this into 4 control volumes centred on the nodes of the 
original mesh). This level of mesh resolution is considered to be the minimum that is acceptable for 
this type of core representation. Increasing the number of cells across the thickness would increase 
the level of resolution of the inter-wrapper flow and could be considered if this model is developed 
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in the future (although it may make the mesh interface with the pools above and below the core 


harder to implement). 


Core Barrel, Reactor Vessel and Guard Vessel 


Since the thermal conductivity of the main structural components is of the same order of magnitude 


as the lead coolant, the temperature gradients across the thickness of the structural components 


are not expected to be significantly different from those within the fluid. To reduce the total cell 


count, the core barrel, RV and GV were meshed using 3 or 4 hexahedral cells through the wall 


thickness (the mesh for the core barrel can be seen in Figure 4.9). Hexahedral cells were used 


to maintain high cell quality and align the cells with the direction of the predominant temperature 


gradient (through the thickness of the wall). 


It is noted that some of the solid components in the LFR have a complex three-dimensional shape, 


such as the support beams that connect the core barrel to the RV. The surfaces of some of the 


solid components are also subject to heat fluxes that vary significantly across the surface. For 


example, the heat flux on the outer surface of the GV is much higher below the water level than 


above the water level. For these components, 3 or 4 hexahedral cells across the thickness of 


the solid component may not be sufficient to accurately compute the local temperature gradients. 


If the CFD model is used in a future analysis to compute the local temperature distribution of 


key structural components, then further refinements of these components may be a necessary 


future development. Further refinement may also be required if the CFD model is used for a future 


analysis where the heat fluxes vary rapidly in time, as this may lead to temperature gradients that 


require increased near-surface mesh resolution in the solids to capture them accurately. 


Material Properties 


Solid Components 


In order to carry out the CFD analysis, the solid components in the reactor (e.g. core barrel, RV 


and GV) require the following material properties to be defined over the range of temperatures 


encountered in the analysis: 


¢ Thermal conductivity. 
* Specific heat capacity. 


* Density. 

Property Correlation Applicable Range 

p 1000(8.0842 — 4.2086 x 10-*7 — 3.8942 x 10°-®T?) 300K < T < 1600K 
on 4186.8(0.1097 + 3.174 x 10-°T) 298K < T < 1670K 
k 100(9.248 x 10-2 + 1.571 x 10-4T) 298K < T < 1670K 


Table 4.4: Material properties of Stainless Steel 316 (SS316) that were used in the CFD 
analysis (Kim, 1975). All quantities are in SI units. 


The majority of solid components in the reactor in this case study are made of stainless steel and 


the material properties are well defined over a large range of temperatures. The material properties 
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for Stainless Steel 316 (SS316) were defined in CFX to vary with temperature and are shown in 
Table 4.4. 


Gaseous Regions 


For gaseous regions (the RV-GV cavity), the following material properties are required over the 
range of temperatures encountered in the analysis: 


¢ Thermal conductivity. 

* Specific heat capacity at constant pressure or volume. 

* Dynamic viscosity. 

¢ An equation of state to calculate the gas density as a function of other state variables. 

¢ Absorption coefficient, scattering coefficient and refractive index, if the gas participates in 
radiation. 


Argon is a gas where reliable material properties are readily available. In this case study, the 
material properties were extracted from the NIST Chemistry WebBook (NIST, 2020) over a wide 
range of temperatures (100°C to 700°C) and were included in the CFD model as an interpolation 
table. 


The ideal gas law was selected for the equation of state to calculate the density of argon. 


In the RV-GV cavity, the majority of heat transfer is expected to occur by radiation and the argon 
will not participate in the radiation. Therefore an absorption coefficient, scattering coefficient and 
refractive index are not required for the simulation. 


Liquid Regions 


For liquid regions (the lead coolant), the following material properties are required over the range 
of temperatures encountered in the analysis: 


¢ Thermal conductivity. 
* Specific heat capacity. 
* Density. 

* Dynamic viscosity. 


For liquid metals, the material properties stated above are often not well defined over the range 
of temperatures encountered in the analysis. This uncertainty in the material properties of liquid 
metals is discussed further in Volume 5. 


Property Correlation Applicable Range 

p 11441 — 1.2795T 600K < T < 1500K 

on 176.2 — 4.923 x 10-7 +. 1.544 x 10-°T* —1.524x 10°T-*. 600K < T < 1500K 

k 9.2+0.011T 600K < T < 1300K 
4.55 x 10-* ettee/7 600K < T < 1473K 


Table 4.5: Material properties of liquid lead that were used in the CFD analysis (NSC, 
2015). All quantities are in SI units. 


Material properties of lead were taken from the recommended correlations in the Handbook on 
Lead-bismuth Eutectic Alloy and Lead Properties (NSC, 2015). These properties are shown in 
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Table 4.5. As noted in the handbook, there is a degree of uncertainty in the material properties of 
liquid lead, with different authors proposing different curve fits. The curve fits adopted in this case 
study are those defined as best estimates in the handbook and are considered adequate for this 
early stage analysis. 


Model Setup 


Turbulence Modelling 


The nature of the case study is such that a RANS approach to turbulence modelling is considered 
appropriate. As outlined in Volume 5 (Liquid Metal Thermal Hydraulics), RANS turbulence models 
have shortcomings when used for the modelling of buoyancy influenced flow and liquid metals, but 
a higher fidelity approach is beyond the resource constraints of this case study. Starting with a 
RANS approach is an appropriate first step even in cases where the intention is to proceed to LES 
modelling in the future. 


The k-w SST model, a standard linear eddy viscosity model, has been chosen to model turbulence 
in this case. While the accuracy of the k-w SST model for natural convection driven flows of liquid 
metal is uncertain, it was found to be the most suitable for buoyancy driven flow from the range of 
models considered by Study A (Liquid Metal CFD Modelling of the TALL-3D Test Facility). It is also 
known to be reasonably stable and robust and is likely to form a good starting point for the early 
stage calculations considered in this case study. 


As discussed in Volume 5 (Liquid Metal Thermal Hydraulics), a simple RANS turbulence model 
may need to be modified when used to simulate liquid metals to avoid excessive turbulent transport 
of thermal energy. In the absence of specialist turbulence models, the simplest modification is to 
change the turbulent Prandtl number, which is an option that is readily available in most CFD codes. 


Increasing the turbulent Prandtl number reduces the turbulent diffusivity of thermal energy relative 
to the turbulent diffusivity of momentum. Ideally this modification would be informed by using ex- 
perimental data from benchmark studies, such as the example considered in Study A (Liquid Metal 
CFD Modelling of the TALL-3D Test Facility). For this case study, specifically relevant information 
is not available, so it is not clear what value of turbulent Prandtl number to use. 


Drawing from previous work (as discussed in Section 3.5), a value between 0.85 (the default in the 
CFD code) and 4 is judged to be most likely to be appropriate for liquid lead (as recommended by 
Roelofs, 2019). In this case study, a value of 2 was chosen for the baseline CFD model. Additional 
simulations were also carried out with turbulent Prandtl numbers of 0.85 and 4, to investigate 
the sensitivity of the model predictions to turbulent Prandtl number. The findings of this study are 
presented in Section 5.4.2. 
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Near Wall Modelling 


Different CFD codes offer their own approaches for modelling the near wall flow. These approaches 
can differ significantly for low Prandtl number fluids (like liquid metals). Therefore, an engineer 
should consult the manual for their CFD code to check the implementation of near wall modelling 
approach that they are adopting. 


In this case study, the ‘Automatic Wall Treatment’ offered by CFX was used to model the near 
wall flow, which is described in the ANSYS CFX Theory Manual (ANSYS, 2020). The Automatic 
Wall Treatment approach is designed to work well across a wide range of y* values, by blending 
together the viscous / conduction sublayer and the logarithmic regions of the momentum and ther- 
mal boundary layers. The final equation for T* is a function of both the dimensionless wall distance 
and the local molecular Prandtl number. This allows the thickness of the conduction sublayer to dif- 
fer from that of the viscous sublayer, as required for the modelling of liquid metals. Furthermore, 
the ‘Automatic Wall Treatment’ can be applied to the near-wall region in both the RV-GV cavity 
(y* < 3) and the lead pools (50 < yt < 100), as the blending depends on the local molecular 
Prandlt number (which is significantly different for argon and lead). 


Buoyancy 


For this CFD model, buoyancy forces are calculated directly from the local density differences using 
the ‘Density Difference’ model in CFX, rather than the Boussinesq approximation. As discussed in 
Volume 3 (Natural Convection and Passive Cooling), the Boussinesq approximation is known to be 
inaccurate when the density changes (relative to the mean density of the fluid) are significant. In 
this case study, the local temperature (and therefore density) variations within the pools and within 
RV-GV cavity could be significant and the ‘Density Difference’ model is more appropriate. 


The effect of buoyancy on turbulence was accounted for by including a source term for the pro- 
duction of turbulent kinetic energy in the transport equation for k, as recommended in the ANSYS 
CFX user manual (ANSYS, 2020). The effect of buoyancy on the specific dissipation rate w is less 
well understood. Therefore, the additional source term was not included in the transport equation 
for w. Sensitivity studies to assess the influence of a range of turbulence modelling choices could 
be carried out in the future, but are not demonstrated in this case study. Examples for different 
applications are included in Case Studies A and B. 


Numerical Settings 


The following numerical settings were adopted for the discretisation and solution of the momentum, 
turbulence and energy equations: 


* Second-order discretisation schemes (termed ‘High Resolution’ in CFX) were used for the 
convection and diffusion terms. Second-order discretisation schemes were used in prefer- 
ence to first-order schemes, as first-order schemes are likely to artificially smear flow struc- 
tures that develop in the lead pools and RV-GV cavity. 


* Second-order backwards differencing (termed ‘Second Order Backwards Euler’ in CFX) was 
used for the temporal derivative in the unsteady simulation. 
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* CFX addresses the pressure-velocity coupling by solving the momentum and continuity 
equations together with a coupled matrix approach. The option to select other pressure- 
velocity coupling schemes (like SIMPLE or PISO) is not available in CFX, so the default 
matrix coupling approach was used. 


The steady-state solver in CFX uses pseudo time-stepping to approach the steady-state 
solution. The pseudo time step size was reduced by a factor of 5 from the default value 
calculated by CFX to obtain a stable solution in the lead pools. This reduction factor was 
obtained by trial and error. 


Solution Procedure 


For natural convection driven flows with relatively weak buoyancy forcing, CFD codes can some- 
times achieve a well converged ‘steady-state’ solution. However, when the buoyancy forcing is 
strong, it may not be possible to achieve a ‘steady-state’ solution, as the unsteady flow features 
are significant and can lead to instabilities in the solution. 


For new analyses, such as the natural convection driven flow in this case study, it is difficult to 
predict whether a ‘steady-state’ solution can be achieved with the CFD solver. Furthermore, the 
flow may be unstable and challenging to converge, unless a good initial condition is provided for 
the flow field. 


The analysis in this case study is challenging to converge and attempting to solve the complete 
model from uniform initial conditions was not successful (i.e. the solution diverges), so the solution 
was developed in stages. Each stage is intended to improve the initial condition that is provided to 
the CFD solver by sequentially increasing the complexity of the model. The requirement to build 
up the complexity of the solution in this way is common to many models and is extremely case 
specific (i.e. the method below has been developed for this particular model)®. A degree of trial 
and error is often required for a completely new model, but the principle of starting with a reduced 
(or less coupled) set of physics and potentially first-order discretisation schemes and then adding 
in additional complexity one step at a time, is common to all. 


For this model, the solution was developed in the following stages: 


1. Start by applying an artificial momentum source to drive the flow around the primary circuit 
loop in the direction expected. The momentum source was applied inside the pump housing 
at the likely location of the impeller, with a magnitude calculated from the natural circulation 
rate predicted by the system code analysis. For this initial stage, the influence of buoyancy 
was turned off, by setting the lead density to a constant value (rather than being temperature 
dependent) equal to the density of lead at the bulk (volume-averaged) lead temperature 
from the results of the system code analysis. This stage in the analysis was found to be 
significantly easier to converge than the final natural convection driven flow field and was 
also used to carry out preliminary checks of the model setup and mesh quality. 


2. Introduce buoyancy forces in the lead by enabling temperature dependent density. This stage 


8 Building up the solution in this way requires assumptions to be made regarding the likely final flow patterns. Where these 
patterns are uncertain, it is possible to ‘lead’ the solution in an inappropriate direction. Care should be taken, especially if it 
is thought possible that more than one stable flow solution exists. 
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enables the temperature and density fields to develop in the lead, while minimising the insta- 
bilities developed by enabling buoyancy forcing. 


3. Gradually remove the artificial momentum source and allow the natural convection driven 
flow field to develop from the forced flow solution. This is the stage that is most likely to 
diverge, and requires careful monitoring of the flow field and equation residuals. 


4. Once the artificial momentum source has been removed completely, examine the ‘steady- 
state’ flow field calculated by the CFD solver (noting that this solution is not a true steady- 
state solution). This flow field may be poorly converged, as the CFD solver uses large ef- 
fective time steps to approach a steady-state solution. If necessary, continue to analyse the 
solution to reach convergence. 


5. In this case study, unsteady flow features were found to be significant in both the lead pools 
and RV-GV cavity. Therefore a well converged flow field was not possible using an entirely 
steady-state approach and a short transient simulation was carried out. The steady-state 
solution was used as an initial condition for the transient to reduce the simulation time. 


The stepwise solution procedure also allows some types of error or inadequate assumptions to be 
identified before too much computational time has been spent and self-checking of the solution is 
recommended at each stage. For example, the specification of the pressure losses coefficients in 
the porous media is most easily checked at the ‘forced convection’ stage of the solution. 


For this CFD analysis there was some concern regarding the closed nature of the solution domain 
(i.e. no inlets of outlets)’. It was anticipated that this could lead to solution difficulties, because the 
volume of the domain is fixed but the lead was expected to have a variable density, so not to be 
restricted to using the Boussinesq approximation. When lead is heated, its density reduces, which 
will very slightly reduce the mass of lead in the CFD model, because the volume is fixed by the 
closed system. It is noted that a very small reduction in the total mass of lead would not be expected 
to make any noticeable difference to the model solution, but could cause convergence difficulties. 
In anticipation of these difficulties, a small aperture was initially included at the top of the lead pool, 
which would allow a small flow of lead to pass in and out of the domain and mitigate the change 
in mass. However, when this aperture was included, the velocities in the vicinity of the aperture 
were found to be unrealistically high. Therefore, the aperture was left closed (as part of the no- 
shear wall) in the final model. The staged development of the solution, in addition to appropriate 
initialisation of the lead temperature and reductions to the artificial time step, allowed solution 
stability to be maintained. The mass imbalance in the lead regions was carefully monitored along 
with the equation residuals during convergence, to check that it did not indicate large discrepancies. 


Case List 


The aim of this case study was to demonstrate the level of analysis that can be achieved within a 
time frame that is representative of early stage reactor development (say 10 weeks). As the CFD 
model is reasonably complex to setup, run and post-process, only a limited number of cases are 
possible in this time period. Therefore, cases were carefully selected, with a focus on the calcula- 
tions that are likely to provide the most useful information to an early stage reactor development 


This issue is described and discussed in more detail in Case Study B: Fuel Assembly CFD and UQ for a Molten Salt 
Reactor, including the use of a reference pressure in a closed domain. 
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process. Table 4.6 summarises the cases that were run with the CFD model. 


Case Simulation Type Pre Approx. Calculation Time 
1 (Baseline) Steady-state 2.0 4 days 

2 Steady-state 0.85 12 hours 

3 Steady-state 4.0 12 hours 

4 Transient 2.0 2 weeks 


Table 4.6: A summary of the cases that were simulated with the CFD model. The ap- 
proximate calculation time gives the length of time for the CFD solver to complete the 
calculation on 40 parallel processors. 


The transient (unsteady) simulation was found to require a significant duration of simulated time 
to resolve the unsteady flow features. This is likely to be due to the size of the domain, the low 
flow velocities and the inertia of the lead (significantly higher than the inertia of traditional coolants 
because of its high density). Therefore, only a single transient simulation was carried out, to in- 
vestigate the differences between a steady and an unsteady approach. For the comparison of the 
different turbulent Prandtl numbers, the steady-state solutions were used. 


Temporal Resolution 


For the transient simulation, a fixed time step of 0.01 s was adopted. This time step was found to be 
sufficient to achieve a Courant-Friedrichs Lewy (CFL) number less than 1 throughout the majority 
of the domain. As the flow is driven by natural convection, a CFL number less than 1 is desirable 
to resolve the unsteady flow patterns. Whilst a larger time step could have been adopted, this may 
lead to a reduction in solution accuracy, as larger CFL numbers indicate that the flow is travelling 
further (through more than 1 cell) each time step. Additionally a higher CFL number can introduce 
instability in the solution and may lead to divergence of the solver. 


The transient simulation was run for a total of 100s of simulated time (10,000 time steps). The 
adequacy of this length of simulation is discussed in Section 5.3. 


Monitoring Convergence 


To track the progress of a CFD model and to assess the level of convergence of the solution, 
it is important to record some information throughout the solution procedure. Most CFD codes 
recommend monitoring the equation residuals, along with integral parameters such as the total 
heat flux and mass flow rates. The global maximum and minimum values of velocity, temperature 
and pressure can also be useful to check that the solution is within sensible bounds. Giving some 
thought to the results of importance and the expected behaviour of the solution is also useful at 
this stage. 


Particular care should be taken where the timescales of the system response vary considerably 
between the different areas of the modelling domain. For example, the temperature of a structural 
component being heated slowly by an adjacent flow might appear to have reached convergence 
when examining the energy equation residuals (i.e. the temperature is only changing by a very 
small amount between one iteration and the next). However, specific monitoring of that component 
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may reveal that the temperature is still changing progressively and has not actually reached a 
steady value. 


In this case study, the following quantities were used to monitor the convergence of the steady-state 
calculations: 


¢ Point monitors for velocity and temperature in the argon-filled RV-GV cavity and lead pools. 

¢ Integral monitors for the total and radiative heat flux across the argon-filled cavity. 

« Amass flow rate monitor for the circulation rate of lead passing between the lead pools. 

¢ Monitors for the volume-average temperature of the large solid components (core barrel, RV 
and GV). 

« Maximum and Root Mean Square (RMS) of the transport equation residuals. 


With unsteady flow features developing in the lead pools, it was not possible to achieve a fully 
converged ‘steady-state’ solution with constant values for the equation residuals and flow field 
monitors. All residuals and point monitors were found to oscillate around a mean value. Further- 
more, the maximum and RMS residual for the momentum equation were found to be several orders 
of magnitude higher than those of the other transport equations. These observations indicate that 
a transient solution is required to achieve complete convergence of the unsteady flow features. 


For the transient simulation, the monitors from the ‘steady-state’ calculations were retained to mon- 
itor convergence. In addition to these monitors, the flow field was exported every 5 seconds of 
simulated time (one twentieth of the total simulation time) to enable a more detailed examination 
of the progress of the solution. 
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The post-processing of a complex CFD model is a significant undertaking and needs to include 
considerable scrutiny of the model predictions regarding how closely it is believed that they rep- 
resent reality, as well as consideration of how to clearly present the findings to interested parties. 
In this case study scenario, it is likely that this CFD model will be the first of many to contribute 
to the reactor development. Therefore, learning from the model and making recommendations for 
improvement to future models is also of considerable value. 


The current CFD model is extensive and could be processed to give predictions of many different 
aspects of the thermal hydraulic performance of the reactor. To maintain the readability of this 
document, only a selection of the available results are presented and discussed as examples. 
Other predictions could be extracted from the model and interrogated in a similar manner. 


The primary results of importance from the CFD model were described in Section 3.1 as follows: 


* Primary coolant temperature distribution within and between the pools. 
¢ Identification of areas of stagnation and/or stratification within the primary circuit. 
¢ Primary coolant natural circulation flow rate. 


¢ Heat transfer across the cavity between the RV and the GV. 


The first three of these results are purely concerned with the flow and temperature distribution 
within the primary circuit lead, while the fourth (heat transfer across the cavity between the RV and 
the GV) requires a more global consideration of the heat transfer in different regions of the reactor. 


The baseline model predictions of the flow and temperature distribution in the primary circuit lead 
are presented in Section 5.1. The predictions of heat transfer from the primary circuit lead to the 
PHRS heat sink are explored in Section 5.2. 


It is not possible to quantitatively assess the adequacy of the model with a high level of confidence 
without the availability of validation data. However, the potential impact of the modelling method 
on the predictions has been considered to the degree possible within the time and budgetary 
constraints of the case study. The following aspects of the modelling method have therefore been 
considered alongside the interpretation of the steady-state results: 


* The adequacy of the CFD mesh in key areas. 
* The potential significance of geometrical simplifications. 
* The potential significance of the circumferentially periodic boundary conditions. 


To improve confidence in the model and to focus future modelling improvements, further interroga- 
tion of the CFD results and sensitivity studies are presented in Sections 5.3 to 5.4: 
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« Acomparison between the steady and unsteady solution results. 

¢ Acomparison of the CFD results with the existing system code analyses. 

¢ An initial investigation into the applicability of the chosen turbulence modelling approach to 
the primary circuit coolant. 


The Primary Circuit Lead 


Overall Flow and Temperature Distribution 


In the early stages of interpreting the model, it is useful to observe the overall flow and temperature 
distribution throughout the reactor first, to gain general understanding of the model predictions 
before examining specific areas or phenomena in detail. 


Figures 5.1, 5.2 and 5.3 show contour plots of the temperature distribution on vertical and hori- 
zontal slices through the model domain, which are useful for making qualitative observations and 
improving understanding of the level of complexity of the flow and temperature fields. Figure 5.1 
also includes a vector plot illustrating the velocity magnitude and direction of the flow on a vertical 
plane passing through the centre of the heat exchanger. 


The following features of the model predications suggest that it is performing well: 


¢ The flow in the natural circulation loop is continuous. 
¢ The temperature and velocity fields are plausible and vary smoothly. 


¢ There is no flow into or through locations that were intended to be solid surfaces and there 
are no ‘walls’ present that unintentionally block the flow. 


¢ The lead moves upwards when it is heated and moves downwards when it is cooled, indicat- 
ing that gravity is turned on and is acting in the correction direction. 


¢ The flow into and out of the porous media is in a sensible direction (e.g. upwards through the 
core and downwards through the pump). 


* The maximum and minimum temperature and velocity are within a sensible range of the 
mean values. 


Additionally, the heat flux into the primary circuit from the fuel balances (to within 0.3%) the total 
heat flux out of the primary circuit through the RV wall and top surface of the pools. 


The highest predicted temperatures occur above the core in the hot pool and the lowest tempera- 
tures occur beneath the core in the lower cold pool. With the exception of the cold region beneath 
the core, the temperature differences within the pools are around 20°C. The temperature differ- 
ences within the primary circuit are significantly lower than the temperature difference between the 
lead and the PHRS heat sink (the water pool), which are of an order of magnitude higher (hundreds 
of degrees). The horizontal slices in Figure 5.3 indicate that the temperature difference between 
the hot pool and lower cold pool is approximately 20°C, while the temperature difference between 
the hot pool and the upper cold pool is approximately 10 °C. 


It is noted that the temperature variation within each pool is predicted to be similar to the tempera- 
ture differences between pools. This implies that for a hand calculation or system code analysis of 
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Figure 5.1: Temperature contours and velocity vectors on a vertical slice through the 
centre of the heat exchanger. 
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Figure 5.2: Contours of temperature on a vertical slice through the centre of the pump. 
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Figure 5.3: Contours of temperature on two horizontal slices at (a) z = 2m and (b) z = 5m. 
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the lead pools, it may not be reasonable to approximate the lead pools with a single bulk (volume- 
averaged) temperature. 


The observations of the temperature variation in the lead pools are in line with expectations. The 
core is the only source of thermal energy and the majority of this energy is transferred out of reactor 
primary circuit through the RV wall. An initial look at the overall model predictions suggests that 
they are realistic. This does not mean the predictions are completely accurate, but does improve 
confidence that there are no gross errors in the model setup (this is an important stage of self- 
checking the model). 


The velocity vectors in Figure 5.1 are useful for understanding the overall behaviour of the flow 
in the lead pools. In the hot pool, the flow is driven upwards by buoyancy and is fastest near the 
centreline of the reactor. In the lower cold pool, the flow moves downwards along the RV wall and 
then up into the core. Away from the walls in the large open areas of the lead pools, Figure 5.1 
shows that the flow moves very slowly. 


These observations of the velocity vectors in the lead pools are in line with expectations. The flow 
in the lead pools is driven by natural convection and the temperature differences within the bulk 
of the lead pools are low. It is therefore expected that the flow in the bulk of the lead pools will be 
relatively slow moving. In the same manner as the temperature field, these reasonableness checks 
do not mean that the model predictions are accurate but do improve confidence that there are no 
gross errors in the model. 


(a) Steady-State (b) Transient (time-averaged) 


Figure 5.4: Contours of temperature on a vertical slice through the centre of the heat 
exchanger, calculated with (a) steady-state and (b) transient analyses. 
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Figures 5.1, 5.2 and 5.3 present predictions from the steady-state analysis. An initial comparison 
between the results of the steady and unsteady predictions (shown in Figure 5.4) shows a high 
level of consistency in the overall temperature distribution. It should be emphasised that the un- 
steady (transient) analysis has been performed to examine local unsteady flow phenomena and 
is not a better representation of the postulated reactor fault transient (i.e. it is still subject to the 
limitations of the ‘snapshot’ approach). However, Figure 5.4 suggests that qualitative results from 
the steady-state analysis are not being significantly affected by local unsteady flow features. Dif- 
ferences between the steady and unsteady predictions are investigated further in Section 5.3. 


The identification of areas of the primary circuit where stratification may occur is a primary aim of 
this model, therefore, it is logical to interrogate this aspect of the results in more detail. The three 
main reactor plena: the hot pool, the upper cold pool and the lower cold pool are examined in turn 
below. Unless otherwise indicated, the results are presented for the baseline case (steady-state 
simulation, with a turbulent Prandtl number Pr; = 2.0). 


Hot Pool 


Figure 5.5 shows velocity vectors and temperature contours on a vertical slice through the hot pool. 
The lead temperature is predicted to be highest near the reactor centreline, where the flow from 
the core forms a buoyant plume. The lead then flows radially outwards along the free surface and 
passes down the core barrel as it cools. No significant regions of stratification or stagnation are 
observed. 


Location of Plane 


(a) Plane Location (b) Velocity Vectors (c) Temperature Contours 


Figure 5.5: Velocity and temperature in the hot pool. 
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Within the hot pool, heat is lost from the free surface by thermal radiation and through the core 
barrel wall by conduction, before the lead leaves the hot pool through the heat exchanger. A small 
region of slow moving cooler flow is observed at the bottom outer edge of the hot pool. However, 
this region of flow is unlikely to stagnate, as it is heated by the central plume on the inner side of 
this region and cooled by the core barrel on the outer side of this region. 


The buoyant hot plume occurs around the central axis of the 60° periodic section and may be 
overly constrained by the periodic boundaries in the CFD model. It is also likely that the flow 
distribution (and therefore the temperature distribution) at the exit from the core is not fully repre- 
sentative, due to the lack of information regarding the core inlet and outlet geometry. However, the 
key flow features observed in the hot pool (the presence of one or more buoyant plumes and no 
significant regions of stratification at the bottom of the pool) are still judged to be realistic, despite 
the limitations in the boundary conditions. 


It is also observed in Figure 5.5 that the heat loss from the upper surface of the hot pool appears 
to be promoting mixing near the top of the hot pool. This boundary condition is currently an area 
of uncertainty in the model inputs and therefore the qualitative results that are generated by this 
boundary condition (increased mixing near the top surface) should be interpreted with caution. 
This is one area of the model where a sensitivity study could be used to investigate the potential 
influence of the thermal boundary condition on the flow and temperature distribution in the hot pool. 


Upper Cold Pool 


Figures 5.6 and 5.7 show velocity vectors and temperature contours on two vertical slices through 
the upper cold pool. Figure 5.6 is located at the outlet of the heat exchanger, while Figure 5.7 is 
located between the pump inlet and the heat exchanger outlet (through the bulk of the upper cold 
pool). 


Figure 5.6 is useful for understanding the flow in the upper cold pool because it shows the be- 
haviour of the flow immediately after it leaves the heat exchanger and enters the upper cold pool. 
The flow leaving the heat exchanger is constrained by the RV wall, so although details of the jets 
of fluid leaving the heat exchanger are not resolved by the CFD model, the bulk flow patterns are 
unlikely to be significantly affected by this omission. The heat exchanger is not operational during 
this fault and the flow from the heat exchanger is hotter than the bulk temperature of the upper cold 
pool and is driven upwards by buoyancy as soon as it enters the upper cold pool. 


Figure 5.7 shows that the flow patterns in the pool are complex, with regions of slow recirculation 
throughout the majority of the pool. The figure highlights the influence of the flow pattern on the 
temperature distribution. For example, a circulating region of warmer fluid at top of the pool next to 
the RV, has a visible effect on the temperature contours in this area. These recirculating regions 
are likely to be unsteady and require a transient simulation (and possibly the use of a higher 
fidelity technique, such as LES) to resolve. However, full resolution of these flow features may 
not be necessary for the purposes of this investigation, if the overall mixing and heat transfer are 
reasonably well represented. It is noted that the observation regarding the potential sensitivity of 
solution to the upper surface boundary conditions made for the hot pool is also relevant for the 
upper cold pool and the mixing near the surface may be over-predicted by the current model. 


The temperature contours show that there is a region towards the bottom of the pool that has the 
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Figure 5.6: Velocity and temperature in the upper cold pool at the heat exchanger outlet. 
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Figure 5.7: Velocity and temperature in the upper cold pool. 
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potential to form a stratified temperature distribution. However, as the pump inlet is located at the 
very bottom of the upper cold pool, the cooler fluid is drawn towards the pumps by the overall 
natural circulation flow, so complete stagnation of the flow (or the generation of large temperature 
differences within the pool) is judged to be unlikely. 
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Figure 5.8: Streamlines in the upper cold pool, initiated from the outlet of the heat exchanger. 


To help visualise the flow patterns, Figure 5.8 shows streamlines in the upper cold pool, which were 
initiated at the heat exchanger outlet. Figure 5.8 shows that as the flow moves upwards along the 
wall of the reactor vessel, it cools and moves circumferentially outwards. The flow distribution in the 
bulk of the upper cold pool is complex, as the flow is heated slightly by the wall of the core barrel, 
cooled by the wall of the RV and is drawn towards the pump at the bottom of the upper cold pool. 
The streamlines can be seen crossing back and forth across the periodic boundaries in a number 
of places. 


It is considered unlikely that the resolution of this initial CFD model is sufficient to accurately predict 
the flow patterns or the exact degree of stratification within the upper cold pool. As the flow in 
the upper cold pool is unsteady, a sufficiently long time-averaged transient calculation would be 
required to give detailed predictions of the flow patterns and stratification. The steady-state solution 
used here (and the duration of the transient calculation) are both not yet sufficient to provide these 
detailed predictions. However, the model does illustrate how the overall design of the upper cold 
pool in combination with the PHRS makes the formation of significant areas of stagnated flow 
unlikely under these conditions. 


Figure 5.9 shows velocity vectors on a vertical plane along the centreline of the heat exchanger. 
The velocity vectors show that the flow is predominately parallel with the axis of the heat exchanger 
but not uniformly distributed, with higher velocities towards the top. Near the bottom of the heat 
exchanger outlet, the flow appears to stall and a small region of backflow is predicted. However, 
this prediction requires several caveats: 


¢ The representation of the heat exchanger using porous media may be allowing backflow to 
occur in a way that would not be possible in reality. 


¢ The mesh resolution at the bottom of the upper cold pool is not sufficiently fine to predict this 
flow with a high level of confidence. 
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¢ The accuracy of the solution at the interface between the porous media region and fluid 
region is likely to be low due to the numerical discontinuity in the momentum sink between 
the porous region (non-zero) and the fluid region (zero). 
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Figure 5.9: Velocity vectors on a plane through the centreline of the heat exchanger, 
highlighting the region of reversed flow near the bottom of the heat exchanger. 


It is therefore not advisable to draw any firm conclusions regarding the specific flow patterns in 
this region of the reactor using the current CFD model. However, this observation, along with the 
observations of the overall flow and temperature distribution within the hot pool, does suggest 
that the flow through the heat exchanger may not be evenly distributed under natural circulation 
conditions. It is also noted that improvements in the CFD model in this area, both in terms of the 
mesh and the representation of the heat exchanger, could be made to study this flow in future. 


Lower Cold Pool 


Figure 5.10 shows contours of temperature on a vertical plane through the bulk of the lower cold 
pool, away from the pump exit. Beneath the core, the contour lines are approximately horizontal, 
indicating that stratification is predicted. In this case, stratification arises because the heat sink (the 
cold wall of the RV) is below the heat source (the core). 


Figure 5.11 shows velocity vectors in the lower cold pool, to highlight the natural convection driven 
flow patterns. In the same manner as the upper cold pool, the lead in the lower cold pool is cooled 
by the RV wall, increasing its density. The lead flows down the RV wall (due to its negative buoy- 
ancy) before passing up into the core. As the lead turns up into the core, it bypasses the stratified 
region at the very bottom of the RV. Away from the RV wall, in the bulk of the pool, the flow velocity 
is significantly lower than close to the wall. 
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Figure 5.10: Temperature in the lower cold pool. 
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Figure 5.11: Velocity vectors in the lower cold pool, coloured by velocity magnitude. 
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Figure 5.11 suggests that a higher proportion of the flow is predicted to enter the outer fuel as- 
semblies in the core compared with those in the centre (and this has been checked by examination 
of the predicted flow distribution through the core). However, this prediction should be viewed with 
caution. In the real reactor, it is likely that the structures at the bottom of the core will have a sig- 
nificant effect on both the flow distribution entering the core and the occurrence of stratification 
beneath the core. These structures are not included in the current model. However, a porous re- 
gion has been included below the core to allow additional (and/or non-uniform) flow resistance to 
be included in the future. 


The main flow feature of interest in the lower cold pool is the density-driven current close to the RV 
wall. Considering the thermal hydraulic mechanisms in the lower cold pool, this flow feature is likely 
to be present in the real reactor. However, the resolution of this flow feature in the CFD model is 
likely to be affected by both the mesh refinement and the wall modelling assumptions in the vicinity 
of the wall. 


To assess the resolution of the computational mesh and the applicability of wall modelling assump- 
tions in the vicinity of the RV wall, Figure 5.12 shows the velocity and temperature profiles on 
several lines normal to the RV wall. The lines have been defined using a polar coordinate system, 
with 6 = 0° taken as a horizontal line from the top of the core. 


Figure 5.12 also shows the velocity and temperature profiles normal to the RV wall expressed in 
dimensionless units (commonly used in the consideration of near wall modelling with CFD). The 


dimensionless wall distance y is: 
+ pu;y 
y — 


pL 


where 9 is the lead density, yz is the lead dynamic viscosity, u,; = ,/Tw/p is the friction velocity and 
T,, is the wall shear stress. The dimensionless velocity U* is: 


Ut = U/u, 


The dimensionless temperature T* is: 


where Ty is the local wall temperature, 6; = qw/pCpu; is the friction temperature, qw is the wall 
heat flux, p is the density, c, is the specific heat capacity. 


Figure 5.12 shows that the FCH of 1.5mm in the lower cold pool is sufficient to achieve y* in the 
range of 80 to 100. For simple boundary layers of fluids with Pr ~ 1, such as forced convection 
of air over a flat plate, there are various informal guidelines regarding the recommended values of 
yt and the mesh spacing in CFD models. However, as the nature of the flow under investigation 
here departs significantly from these conditions, engineering judgement will be used to assess the 
mesh resolution in the CFD model. 


Figure 5.12 shows that the temperature increases monotonically with increasing distance from the 
wall. Within ~ 0.1 m of the wall, the lead temperature reaches a temperature representative of the 
bulk pool temperature. It can be estimated that there are 10 to 15 cells between this location and the 
wall. From this observation it appears that the mesh resolution is sufficient to capture the thermal 
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Figure 5.12: Velocity magnitude (left) and temperature (right) along lines from the top 
centre of the core to wall of the reactor vessel in the lower cold pool. The profiles are 
plotted in dimensioned units (top) and dimensionless units (bottom). 
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profile. In addition to providing at least 10 cells across the thickness of the thermal boundary layer, 
it is desirable to size the first cell such that the ANSYS CFX Automatic Wall Treatment assumes an 
approximately linear temperature profile (i.e. well within the conduction sublayer)' between the cell 
and the wall. Consideration of the blending function for the Automatic Wall Treatment (as defined 
in ANSYS, 2020), suggests that y* < 90 should be sufficient for the wall function to achieve a 90% 
linear temperature profile across the first cell. 


In contrast to the temperature, the velocity magnitude does not increase monotonically, but reaches 
a peak within ~ 0.01 m of the wall. The magnitude of the peak velocity is ~ 0.15 m/s, which is 
approximately 5 times the velocity magnitude far away from the wall (~ 0.03 m/s). The peak velocity 
occurs at y* ~ 300 to 500 and there are 3 to 4 cells between the peak velocity and the wall. 
Examination of the velocity profile suggests that this is not enough and it is judged that around 8 to 
10 cells between the location of peak velocity and the RV wall would be more appropriate for good 
resolution of the velocity profile. 


As expected for a fluid with a low molecular Pr, the thermal boundary layer appears considerably 
thicker than the momentum boundary layer. However, it is not possible to define exact momentum 
and thermal boundary layer thicknesses, as it is difficult to define a freestream velocity and tem- 
perature for the lead pools. Despite this limitation, the analysis presented in this section suggests 
that the thermal boundary layer is likely to be sufficiently resolved with the current mesh, while 
the momentum boundary layer would benefit from mesh refinement. The low molecular Prandtl 
number results in an increase in the thickness of the conduction sublayer and therefore a reduction 
in the mesh resolution required for the temperature profile relative to the velocity profile. A mesh 
sensitivity study (potentially using a much smaller sub-model) could be performed to quantify the 
changes in the wall shear stress and wall heat transfer prediction that result from refinement of the 
mesh, to help define a sufficient resolution for future models. 


Figure 5.13 shows contours of temperature and velocity magnitude on a plane passing through the 
pump. The density-driven current can still be seen flowing along the RV wall and it does not appear 
to be significantly disrupted by flow from the outlet of the pump housing. It can be observed that 
higher velocity flow from the pump outlet penetrates a short distance into the lower cold pool. How- 
ever, its warmer temperature (compared with the bulk of the lower cold pool) results in buoyancy 
forces opposing the downwards direction of the flow entering the lower cold pool. The buoyancy 
forces reduce the velocity of the downwards flow entering the pool and then directs it upwards 
towards the top of the pool. 


It is noted in examination of the flow predictions in this pool, that the qualitative temperature distri- 
bution may be significantly altered if the water level in the PHRS heat sink were lower. Considering 
a situation where the water level is below the pump outlet, the density current may not be present 
(or may be much weaker) above the height of the water. Under these conditions it may be possible 
for increased stratification to occur in the top half of the lower cold pool as the warmer lead enters 
and then rises. This could be investigated with this model in the future. 


1 Reasons for this are discussed in Volume 5 
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Figure 5.13: Contours of temperature and velocity on a plane through the pump in the lower cold pool. 
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Primary Coolant Natural Circulation Flow Rate 


Prediction of the overall natural circulation flow rate within the primary circuit lead was a specific 
objective of this analysis and, as such, the flow rate of lead through the pump was monitored 
throughout the solution. However, it was noted in the steady-state analysis that this was an aspect 
of the simulation that did not reach a stable value. 
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Figure 5.14: Circulation rate predicted by (a) steady-state and (b) transient analyses. 


Figure 5.14a shows the CFD predicted circulation rate over a period of 5,000 iterations for the 
baseline steady-state analysis. It can be seen that the flow rate varies between 1,640kg/s and 
1,820 kg/s with no indication that it is converging towards a stable value (despite the residuals in the 
momentum equations reaching stable values well before this point). This implies that the natural 
circulation flow rate is potentially affected by unsteady flow features and so further investigation 
using the transient analysis results is needed. 


When interpreting the CFD results for the natural circulation rate, the mean circulation rate, the 
magnitude of the unsteady variations around the mean and the nature of the variations (regular 
or chaotic) are of interest. To make an assessment of these features of the natural circulation 
rate, Figure 5.14b shows the predicted variation in natural circulation flow rate with time over 100 
seconds of the transient calculation. The transient calculation predicts that the flow rate varies from 
the initial steady-state snapshot solution of + 1750 kg/s to little below 1650 kg/s before rising again 
to + 1715kg/s. 


Both the steady-state and transient solutions in Figure 5.14 indicate that the mean circulation rate 
is likely to be around 1700 kg/s. However, if the magnitude and nature of the unsteady fluctuations 
around the mean circulation rate are to be assessed, the length of the transient calculation needs 
to be significantly longer than the 100s considered in this analysis (limited to align with the time- 
frames, budgets and computing resources assumed by this case study scenario). Considering the 
predictions of overall circulation rate, it would certainly be of benefit to extend this transient in the 
future to see if the mean circulation rate changes significantly from the steady-state solution, and 
to assess the magnitude and nature of the unsteady variation around the mean. 


Due to the strong coupling between the thermal and momentum fields in a naturally circulating 
flow and the complexity of the reactor primary circuit, the natural circulation flow rate is a model 
prediction where it is difficult to have a high level of confidence in the quantitative result. All of 
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the sources of pressure loss in the circuit (from components or flow phenomena) and all of the 
heat transfer mechanisms potentially have an effect on this prediction, so uncertainty in almost any 
aspect of the model will result in uncertainty in the overall flow rate. 


Due to the modelling assumptions described above, The precise cause of the unsteady nature of 
the natural circulation flow rate is also difficult to identify, although potential candidates are explored 
further in Section 5.3. 


Heat Transfer to the PHRS 


In addition to investigating the primary circuit flow and temperature distributions, this CFD model 
aims to give further insight into the operation and performance of the PHRS. For the PHRS to 
remove heat from the reactor primary circuit, the heat must first be transported through the lead to 
the walls of the RV. The heat must then be transported across the RV-GV cavity, through the GV 
walls and finally into the PHRS heat sink. Due to the path followed by the heat transfer, this part of 
the PHRS can be conceptualised as a series of thermal resistances in series. 


For steady-state conditions, with no other sources of decay heat removal, the rate of heat transfer to 
the PHRS heat sink is equal to the rate of heat generation within the reactor core. The model (and 
indeed the real system) under these conditions will heat up until the temperature difference across 
each component is sufficient to transfer thermal energy at this rate. Since the largest temperature 
difference will occur across the component with the highest thermal resistance, this component 
has the greatest influence on the primary circuit lead (and therefore the reactor fuel) temperatures 
during the fault. 
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Figure 5.15: Approximate temperatures of the key components in the PHRS and the 
conceptual model for this system as a set of thermal resistors in series. 
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Figure 5.15 shows approximate temperatures for key components in the heat transfer path from 
the primary circuit to the PHRS heat sink at two locations; one above the heat sink water surface 
and one below. The temperature gradient between each location in Figure 5.15 is indicative of 
the thermal resistance provided by that component. It is interesting to observe that the effective 
thermal resistance between the primary circuit lead and the RV wall is reasonably low, both above 
and below the water level. The effective thermal resistance between the GV and the PHRS heat 
sink is negligible below the water level but represents the highest resistance above the water level. 
The heat transfer across the RV-GV cavity has a significant influence on the PHRS performance, 
regardless of the water level. However, at the very high temperatures seen at this stage in the 
fault transient, the thermal radiation heat transfer across this cavity is high enough to result in an 
effective thermal conductivity similar to the steel vessels. 


Two model predictions of interest are considered in more detail below: 


* Heat transfer across the RV-GV cavity. 
¢ The RV temperature distribution. 


Results relating to other components could be examined in a similar way, although the heat transfer 
from the GV to the PHRS heat sink is specified by boundary conditions (i.e. it is largely a function 
of the model inputs) and so a detailed interrogation of this model prediction is not likely to yield any 


further insight. 
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Figure 5.16: Vertical profiles of radiative heat flux and total heat flux (radiative and con- 
vective) on the surface of the RV in contact with the argon. The heat fluxes are circum- 
ferential averages. 


Figure 5.16 shows the variation of the radiative and total heat flux (radiative and convective) along 
the height of the RV. The difference between the radiative heat flux and the total heat flux is less 
than 10%. This indicates that: 


¢ The majority of the heat transfer across the cavity occurs by radiation. 


2 Although there is not a high level of confidence in the heat transfer boundary condition specified between the GV and the 
PHRS heat sink above the water level. 
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¢ The emissivity of the outer surface of the RV and the inner surface of the GV may have 
a significant effect on the thermal resistance of the cavity. Care should be taken to ensure 
that the emissivity of these surfaces are known, taking into account their likely ‘in service’ 
condition. 


As the proportion of the total heat transfer that can be attributed to convection is small, it can be 
concluded that accurately predicting the precise convection pattern within the RV-GV cavity is not 
likely to significantly affect the results of interest in this model. However, if the emissivity of the cavity 
surfaces were reduced, then a greater proportion of the heat flux would be transported across the 
cavity by convection. Under these conditions, it would be more important to accurately resolve the 
natural convection driven flow within the cavity with the CFD model. Convection is also likely to be 
a proportionally more important heat transfer mechanism when the absolute temperatures of the 
surfaces are lower (e.g. during normal operation). 
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Figure 5.17: Vertical profiles of wall temperature on the outer surface and the inner 
surface of the RV and the GV. The temperatures are circumferential averages. 


Figure 5.17 shows the variation in the (circumferential average) wall temperature with height on 
the inner and outer surfaces of the RV and GV. These profiles are useful for understanding the 
variation in heat transfer across the cavity with height (Figure 5.16). The temperature of the outer 
surface of the GV is significantly higher above the water level (z > 5 m) than below the water level. 
The temperature is higher above the water level because heat transfer from the outer surface of the 
GV is by natural convection and thermal radiation, which results in significantly less effective heat 
transfer than that provided by the nucleate boiling below the water level. As a result of the increased 
temperature of the outer surface of the GV above the water level, the temperature difference across 
the RV-GV cavity is reduced. The reduced temperature difference reduces the heat flux across the 
top region of the cavity and the strength of the natural convection driven flow in this region of the 
cavity, as shown in Figure 5.18. 


Over the duration of the reactor decay transient, the water level reduces. It is expected that at a 
slightly earlier or later time snapshot, the temperature of the RV and GV walls will have a similar 
characteristic shape to Figure 5.17, with the GV wall temperatures being significantly higher above 
the water level. However, as the water level drops below the core, it is possible that the vessel 
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temperatures could drive substantially different primary circuit flow patterns. These primary circuit 
flow patterns could be investigated by carrying out additional analyses with the CFD model in the 
future. 


Lead --- g=5m 


Guard Vessel 


Reactor Vessel 


=a 0.2 0.4 0.6 0.8 1 


x/L(—) 


Figure 5.18: Horizontal profiles of vertical velocity (W) in the argon-filled cavity. x/L is 
the fractional distance across the cavity. x/L = 0 represents the RV wall and x/L = 1 
represents the GV wall. 


An assessment of the mesh (not shown here) showed that y* ~ 2 on both the RV and GV walls. 
A y* ~ 2 indicates that the cell in contact with the wall is placed well within the viscous and 
conduction sublayers of the momentum and thermal boundary layers in the argon. Hence, the 
momentum and thermal boundary layers are both deemed to be well resolved with this level of 
mesh resolution. Moving away from the first cell into the bulk of the cavity, the velocity increases at 
first and then decreases after reaching a peak at y* ~ 10. The overall mesh resolution provides 
a smooth velocity profile prediction with enough cells between the peak velocity and the wall for 
a reasonably accurate prediction of the wall shear stress and wall heat flux. For the current CFD 
model, where the majority of the heat is transferred across the cavity by radiation, further mesh 
refinement is not deemed to be necessary. 


In the majority of the RV-GV cavity, the flow patterns are typical of natural convection driven flow 
between two vertical plates. However, towards the bottom, the cavity is heated from above by the 
lead and cooled from below by the PHRS. This heating and cooling configuration has the potential 
for stratification. To investigate this, Figure 5.19 shows velocity vectors and temperature contours 
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towards the bottom of the RV-GV cavity. 
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Figure 5.19: Vectors of velocity and contours of temperature at the bottom of the RV-GV 
cavity beneath the reactor core. 


The temperature contours in Figure 5.19 show a region of stratification near the reactor centreline. 
This is supported by the velocity vector plot, which shows that the velocity in this region is low 
compared to the natural convection driven flow higher up in the cavity. This stratification configu- 
ration is stable and is not seen to be disturbed significantly by the natural convection driven flow 
higher up in the cavity, both in the steady and unsteady solutions. In reality, it is unlikely that the 
region of stratified flow at the bottom of the RV-GV cavity is of primary concern. The majority of the 
heat is transferred across the cavity by radiation, therefore, neither the overall heat transfer nor the 
temperature of the solid components (RV and GV) will be significantly affected by the temperature 
of the argon gas. For this analysis, it is sufficient to observe the slow moving flow at the bottom of 
the cavity and note this as a feature of the current PHRS design. 


Temperature of the Reactor Vessel Wall 


The temperature distribution on the RV wall was identified as a result of future interest for this 
CFD model and large temperature gradients could have implications for the structural integrity of 
the RV. While the model is not sufficiently developed to predict the detailed (and localised) vessel 
temperatures with a high level of confidence, the overall RV temperatures directly affect the heat 
transfer across the RV-GV gap and therefore the performance of the PHRS. 


Figure 5.20 shows the temperature of the inner and outer surfaces of the RV wall. Over the majority 
of the RV inner surface (in contact with the lead), the temperature varies between 680°C and 
740°C. Meanwhile, over the majority of the outer surface of the RV (in contact with the argon), 
the temperature varies between 600 °C and 660°C. The temperature difference between the inner 
and outer surfaces of the RV is consistent with the temperature difference required to transport 
the decay heat load across the thickness of the RV (8cm thick stainless steel). The temperature 
difference between the inner and outer surfaces is also in line with the system code predictions, 
which are presented later in Section 5.4.1. 
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Figure 5.20: Temperatures of the inner and outer surface of the RV wall. 


Both the inner and outer surfaces of the RV show a significant drop in temperature at the bottom 
of the reactor, as a result of the thermal stratification that occurs in the lead pool beneath the core. 
However, as discussed in Section 5.1.4, while stratification beneath the core is likely to occur in 
reality, it is also likely to be affected by geometry that is not currently included in the model. Hence, 
the temperature gradient experienced at the bottom of the RV is likely to change in future evolutions 
of the model. 


The influence of the water level can be seen clearly on the outer surface of the RV, but is barely 
noticeable on the inner surface. The temperature distribution on the inner surface of the RV is 
strongly influenced by the adjacent lead, due the the high convective heat transfer coefficient on 
this surface. However, the temperature distribution on the outer surface of the RV is primarily af- 
fected by the temperature distribution on the GV wall, as the outer surface of the RV is cooled by 
radiative heat transfer across the cavity. This observation suggests that three-dimensional thermal 
conduction within the RV is an important mechanism in the prediction of the temperature distribu- 
tion throughout this component. 


The solid components in the CFD model (including the RV and the GV walls) are currently meshed 
with 3 cells across the thickness. This is significantly coarser than the mesh resolution that is 
adopted for the fluid in the lead pools and the RV-GV cavity. The current level of resolution is 
reasonable for an initial model but, given the significant three-dimensional conduction within the RV, 
increased refinement of the mesh through the thickness of the RV is likely to be necessary for future 
analyses. Refinement of the mesh perpendicular to the thickness should also be considered close 
to regions of high vertical temperature gradient. If a similar model were to be used for a transient 
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analysis where the RV was expected to be subject to fluid temperatures that vary significantly in 
time, then the mesh would need to be refined further, as noted in Volume 2 and demonstrated in 
Study D. 


Comparison of Steady and Unsteady Solutions 


The flow within the primary circuit is driven by natural convection so it was expected from the outset 
that there would be areas of unsteadiness due to the presence of inherently unsteady thermal 
hydraulic phenomena. It was therefore anticipated at the planning stage that it may not be possible 
to achieve a converged steady-state solution and that a transient analysis (albeit one with constant 
boundary conditions) should also be performed. 


It is reiterated that this ‘transient’ flow solution is not intended as a subsection of the reactor fault 
transient. The boundary conditions and the starting condition for this analysis both result in the 
solution still being subject to the limitations of a ‘snapshot’ approach. The purpose in this case is 
simply to observe the level of unsteadiness in the flow field and determine if the convergence of 
the momentum equation can be improved if this unsteadiness is properly resolved as a transient 
solution. 


Unsteady fluctuations within the bulk of the lead pools are not likely to be of particular interest. 
However, where this unsteadiness impacts on the results of importance (or those of future interest), 
it needs to be considered. For example, unsteady temperature fluctuations in fluid adjacent to 
structural components may cause fatigue damage. 


In this section, three areas where unsteady flow features would be expected to occur are explored 
and their impact on the results of importance are discussed: 


* Below the core. 
¢ Above the core. 
¢ Upper cold pool. 


Below the Core 


Figure 5.21 shows contours of vertical velocity on a horizontal plane through the lower cold pool, at 
20 second time intervals over the duration of the transient simulation. The contours show significant 
regions of flow unsteadiness, both close to the RV wall and in the bulk of the pool. The contours 
indicate that a point probe placed in the bulk of the lead pool could see variation in both velocity 
magnitude and direction. When compared to the peak velocity near the RV wall (of = -0.16 m/s as 
shown in Figure 5.12), the magnitude of the unsteady fluctuations in this region is large. 


The unsteady fluctuations in the bulk of the pool are driven by the interaction between the following 
forces: 


¢ Buoyancy-driven (upwards) forcing of the hot jet of fluid that enters the pool from the pump 
outlet. 


¢ Buoyancy-driven (upwards) forcing of the flow close to the wall of the core barrel. 


¢ Negative buoyancy-driven (downwards) forcing on the flow close to the wall of the RV. 
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Figure 5.21: Contours of vertical velocity on a horizontal plane through the lower cold 
pool at z = —1.5m. 


* Shear forces between the faster moving jet of fluid that enters the pool and the relatively slow 
moving fluid in the bulk of the pool. 


The magnitude of the unsteady fluctuations that are seen in the bulk of the pool are judged to be 
plausible, given the nature of the buoyancy-driven and shear forces that are acting on the fluid in 
the pool. 


The region of the lower cold pool below the core is likely to be important for both the overall natural 
circulation flow rate and the flow and temperature distribution entering the core. The unsteady 
flow observed in this area appears to be due to the interaction of the density current with the 
Stratification observed at the bottom of the RV. The interaction can be seen in Figure 5.22, which 
shows contours of temperature immediately beneath the core at 20 second time intervals over the 
duration of the transient simulation. 


It is also possible that this interaction contributes to the unsteadiness in the overall natural circula- 
tion flow rate. Although this is an area where the CFD model is missing geometrical features, it is 
likely that there will be some degree of stratification and a density-driven current in the real reactor, 
leading to a degree of unsteadiness in the natural circulation flow rate. 
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Figure 5.22: Contours of temperature on a vertical plane through the bottom of the lower 
cold pool. 


5.3.2 Above the Core 


Another potentially unsteady flow feature that was identified in the steady-state calculation is the 
buoyant plume in the hot pool, which is located near the reactor centreline. To examine this feature, 
Figure 5.23 shows contours of velocity magnitude on a horizontal plane through the hot pool at 20 
second intervals over the duration of the transient. 


The buoyant plume in Figure 5.23 shows some unsteady variation in the location and magnitude 
of the velocity contours over the duration of the transient simulation. However, the changes in the 
location and magnitude of the velocity contours are relatively small. Therefore Figure 5.23 suggests 
that the buoyant plume in the hot pool is likely to be relatively steady and remain in the same 
position. However, as discussed in Section 5.1.2, the buoyant plume may be overly constrained 
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Figure 5.23: Contours of velocity magnitude on a horizontal plane through the hot pool at z = 4m. 


by the periodic boundaries close to the axis of rotation and the plume may ‘wander’ further than 
this simulation predicts. Similar behaviour was observed in the non-stratified jet simulation in Case 
Study A: Liquid Metal CFD Modelling of the TALL-3D Test Facility. 


Furthermore, the above core structure is not included in this CFD model, as details of its design 
were not available at the time of this study. Due to the obstruction caused by the above core 
structure, the flow may form a more complex plume structure than predicted by this CFD model. 
Therefore, if the buoyant plume and flow distribution in the hot pool are of interest, an additional 
analysis using a full 360 ° model (which includes the above core structure) would be required. 


Upper Cold Pool 


The analysis of the flow and temperature in the upper cold pool in Section 5.1.3 suggested that the 
flow field may include unsteady flow features. To examine these unsteady flow features, Figure 5.24 
shows contours of vertical velocity on a horizontal plane through the upper cold pool at 20 second 
time intervals over the duration of the transient. 


As the flow exits the heat exchanger and enters the upper cold pool, it is driven upwards by buoy- 
ancy. The flow then moves circumferentially outwards and is cooled by the RV wall. Figure 5.24 
shows that the upwards (buoyancy-driven) flow at the outlet of the heat exchanger appears to be 
relatively steady. However, the flow towards the middle of the pool, away from the solid surfaces, 
shows significant variation in the location of the areas of high and low vertical velocity over the du- 
ration of the transient. it is noted that the periodic boundaries allow these regions to move across 
them, whereas a symmetry boundary condition would have artificially constrained them 
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Figure 5.24: Contours of vertical velocity on a horizontal plane through the upper cold pool at z = 5m. 


Although the precise variation in flow pattern is not necessarily quantitatively accurate (as the 
duration of the transient simulation has not been not long enough for the initial conditions to have 
fully dissipated), it can be concluded that unsteady flow features are likely to be present in the 
upper cold pool. 


Transient Analysis Duration 


Figure 5.25a shows the variation in the bulk circulation of lead between the pools over the duration 
of the transient simulation, while Figure 5.25b shows point probes of velocity magnitude that were 
placed in the lead pools at the start of the transient simulation. 


The point probes allow the velocity magnitude to be monitored continuously throughout the sim- 
ulation, at the location of the probes. This can be useful because it allows the unsteadiness in 
the flow field to be observed more clearly than the discrete snapshots of the complete flow and 
temperature fields, which were output every 20s of simulated time. However, in a first simulation 
it can be difficult to place the probes in the optimum locations, as the transient behaviour of the 
flow field is unknown before the simulation is carried out. Therefore, in addition to placing probes at 
likely locations of interest, a number of additional point probes were placed at random in the lead 
pools for this first simulation. Having carried out the transient simulation and observed the results, 
the location (and number) of probes can be updated for future transient calculations. For example, 
placing probes at the inlet and outlet of the lead pools would allow an assessment to be made of 
whether the transient oscillations are contained within the lead pools or move between the lead 
pools. 
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Figure 5.25: Transient variation of (a) the lead circulation rate and (b) point velocity 
monitors over the duration of the simulated transient. 


Figure 5.25 shows that the bulk circulation appears to vary with a period of approximately 100 
seconds, with magnitude of around 5% of the mean value. Due to the apparent period being the 
same as the entire length of the simulation, Figure 5.25a suggests that the transient calculation has 
not been run for a sufficient duration to accurately assess the unsteadiness in the bulk circulation 
rate. This long period is likely to arise in part because of the inertia of the lead and the large volume 
of the lead pools. 


In the same manner as the bulk circulation rate, the point velocity monitors in Figure 5.25b have an 
unsteady component with a significant magnitude and a long period (around 50 to 100 seconds). 
The period of the point velocity monitors indicates that the local unsteady flow features within the 
lead pool have a similar period to the bulk circulation rate. Therefore, future unsteady calculations 
are likely to require a longer duration to capture the possible range of both the local flow features 
and bulk circulation rate. 


Despite the short duration of the transient simulation, it has provided several useful outcomes and 
enables: 
¢ Development of some initial understanding of the unsteady flow features within the pools. 


* Confirmation that the major flow features that were identified in the ‘steady-state’ analysis 
are also present in the transient simulation. 


¢ Preparation of a better plan for any future transient analysis, both in terms of numerical 
parameters (like time step and simulation duration) and in the monitoring of flow features of 
importance. 
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Increasing Confidence in the Results 


This analysis has been carried out during the early stages of reactor development so experimental 
measurements are not available to provide specifically relevant validation of the results. Despite 
the lack of validation data, it is valuable to consider other activities to improve confidence in the 
results. 


So far, observations of the CFD results have been used to provide insight into the flow and temper- 
ature fields within the reactor pools and RV-GV cavity. In addition to providing understanding, the 
impact of the periodic boundaries, geometric simplifications and the mesh resolution have been 
considered and the predictions of both steady-state and transient solutions compared. 


For this analysis, the following additional activities have been undertaken to improve confidence in 
the results: 


¢ Independent verification of the model geometry, mesh, boundary conditions, numerical pa- 
rameters for the CFD solver and solution convergence. 


¢ Basic hand calculations to check specific aspects of the CFD model predictions. 


* Acomparison with the results of a system code analysis of the reactor, that had been carried 
out before the CFD analysis. 


¢ A study to investigate the sensitivity of the model predictions to the turbulence modelling 
approach (in this case the turbulent Prandtl number). 


« A study to investigate the influence of the pool top surface thermal boundary condition on 
the model predictions. 


It is noted that the last of these studies was not included in the original analysis plan, but has 
arisen as a result of the observed influence of this boundary condition on the flow in the hot and 
upper cold pools. A number of potential future sensitivity studies have been suggested throughout 
Section 5. However, this one has been selected for immediate consideration as it potentially effects 
judgements regarding the likely level of stratification in the pools (which is a key output from the 
current CFD model). 


The independent verification was carried in line with the plan described in Section 3.9. Basic hand 
calculations are not fully reported here, however, they broadly comprised checking that the: 
* Integral of the heat source within the core matches the intended decay heat. 


* Overall thermal energy balance (i.e. the total heat leaving the model is the same as the heat 
source). 


¢ Heat transfer via different mechanisms (conduction, convection, radiation) through the RV, 
GV and across RV-GV gap are broadly consistent with comparison calculations. 


¢ Pressure drops over key components are consistent with comparison calculations using the 
input loss coefficients. 


The comparison with the system code and investigations of the impact of the turbulent Prandtl 
number and the top surface boundary conditions are presented below. 
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Comparison with System Code Results 


Comparing the results of the CFD analysis with an analysis carried out using a different tool (in this 
case a system code) is useful for several reasons: 


¢ It provides an additional check of the model setup (input parameters and boundary condi- 
tions) since large inconsistencies in the results could be an indication that the CFD model 
has been setup incorrectly. 


* Differences in the integral results from the CFD model and system code results can be used 
to highlight modelling choices, simplifications and assumptions which (whilst not necessarily 
incorrect) require further interrogation. 


¢ Differences in the CFD and system code results can be used to assess the importance and 
impact of three-dimensional flow features. 


However, care should be taken in drawing conclusions about the relative performance and accu- 
racy of the models from a comparison of results from a single case, as differences in predictions 
can originate from a number of different sources. For example, in a research environment such a 
comparison may be made between models that were specifically created for the purpose and are 
therefore set up to represent exactly the same reactor design and fault conditions. Models created 
in an industrial setting are usually aimed at progressing or de-risking aspects of the design, not for 
performing code-to-code comparisons. It is, therefore, not uncommon in a fast moving industrial 
design setting, for differences between models created at different times (often by different engi- 
neers) to be caused by differences in the input data, where each model has attempted to capture 
the latest reasoning at a specific point in the design process. 


When comparing the results of different models, the following features of the models should be 
considered as possible sources of discrepancy: 


¢ Different input data, i.e. are the models actually attempting to represent identical circum- 
stances? 

¢ Differences in user assumptions, set up and analysis procedure. 

* Differences in the underlying principles and approaches within the modelling tool. 


In this investigation, the following parameters are compared between the CFD analysis and the 
system code analysis: 


¢ Volume-averaged temperature of the lead pools. 
¢ Bulk circulation rate of lead between pools. 
* Temperature of the RV and GV walls. 


Temperature of the Lead Pools 


Table 5.1 shows the bulk (volume-averaged) temperature in the lead pools, computed by the CFD 
model and the system code analysis. 
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Region CFD (°C) System code (°C) 
Hot Pool 739 683 
Upper Cold Pool 734 672 
Lower Cold Pool 724 664 
All Pools 730 673 


Table 5.1: A comparison of bulk (volume-averaged) temperature of the lead pools, com- 
puted by the CFD analysis and system code. 


The system code consistently predicts lower bulk temperatures in the lead pools than the CFD 
analysis (~ 60°C lower). The magnitude of this discrepancy is large and needs further investiga- 
tion. A comparison of the input data used by both models revealed several important differences: 


¢ At the time when the system code model was developed, the design used a reactor vessel 
and guard vessel that were both 10% larger than those in the design used for the CFD model. 


* Developments in the design of the PHRS heat sink have resulted in the water level (relative 
to the top of the RV) being lower at the ‘snapshot’ in the reactor decay transient (0.5 days 
after the initial loss of power) in the CFD model than in the system code model. 


¢ Assumptions regarding the variation in core power over the duration of the fault have changed 
such that 0.5 days after the initial loss of power, the heat input in the system code model is 
17% lower than that in the CFD model. 


All of these differences act to increase the bulk lead temperatures in the CFD model relative to 
the system code model. A hand calculation using the ‘thermal resistances in series’ approach 
(described in Section 5.2) can be used to estimate the approximate magnitude of the changes 
in overall predicted lead temperature for each of these differences. In this case it was found that, 
when added together, the differences in these inputs are more than sufficient to account for the 
discrepancy seen in the lead temperature predictions. 


A comparison of the user assumptions and setup of both models found numerous differences, a 
number of which are potentially important to the investigation of the bulk lead temperatures: 


¢ The extent of the two models is substantially different. The system code includes a con- 
crete top shield above the reactor cover gas and a containment above that, in addition to 
components within the PHRS heat sink, none of which are included in the CFD model. 


¢ The level of resolution of the lead pools, RV and GV is considerably lower in the system code 
model than in the CFD model. 


¢ The heat loss from the surface of the lead pools is much greater in the CFD model than in 
the system code. This is largely due to the differences in model extent and the boundary 
conditions applied to this surface in the CFD model. 


¢ The system code model analyses the full fault transient, starting from normal operating con- 
ditions. The CFD model only analyses a ‘snapshot’ of this transient, using a steady-state 
approach. 


The difference in the predicted heat loss from the lead pool surface between these two models acts 
to reduce the bulk lead temperatures in the CFD model compared with the system code model. 
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Accounting for this heat loss in the ‘thermal resistances in series’ hand calculation mentioned 
above, aligns the discrepancy between the bulk lead temperatures predicted by each of the CFD 
and system code models with that which would be expected due to the differences in model inputs. 


It is worth noting that just because the discrepancy in temperature predictions between the two 
models can be explained in this situation by a combination of four differences in input data, does not 
mean that there are no other important differences between the models. As discussed in Volume 4: 
Confidence and Uncertainty, it is quite possible for two important differences to act in different 
directions and cancel each other out, such that they are not apparent in an examination of the 
results of one specific case. However, they may cause larger discrepancies when the models are 
used to analyse a different scenario. 


For example, it was initially considered that the discrepancy in bulk lead temperatures could be 
due to the ‘snapshot in time’ approach taken for the CFD. If the lead in the pools was heating up at 
this stage of the transient then adopting a thermal equilibrium approach (as is the case for the CFD 
model) would result in an over prediction of the bulk lead temperatures. However examination of the 
system code results suggests that this is not the case and the lead pool temperatures predicted by 
the system code are dropping slightly at this time, with the heat flux leaving the model being 1.4% 
higher than the heat entering the model from the core. At this point in the transient, the thermal 
imbalance is small compared with the other differences discussed above. However, this factor is 
likely to be more significant at other points in the transient and could become more important in 
the future when the CFD model is aiming at a greater degree of quantitative accuracy. 


Despite the large difference in overall lead temperature predictions, it is noted from Table 5.1 that 
both the system code and the CFD models predict temperature differences between the lead pools 
of less than 20 °C. When considered in the context of the temperature drop observed between the 
primary circuit lead and the PHRS heat sink, this range is small. 


Bulk Circulation Rate of Lead Between Pools 


Using the results of the transient analysis prediction in Figure 5.25a, the bulk circulation rate of lead 
between pools was calculated to be around 1670 kg/s with the CFD model. In contrast, the system 
code analysis predicted a bulk circulation rate of lead of 2250 kg/s at nominally the same point in 
the fault transient. The significant differences in assumed reactor geometry, core heat source and 
pool surface heat transfer previously discussed are also likely to effect the circulation rate (although 
it much more difficult to perform a simple hand calculation to estimate by how much). 


Differences in the overall lead primary circuit pressure losses will also have an impact on the 
circulation rate. Both the CFD model and the system code analysis assume the same pressure 
loss characteristics for the major components in the primary circuit. However, the CFD model also 
accounts for pressure losses from 3D flow features such as: 


¢ Natural convection driven mixing in the lead pools. 

¢ The heat exchanger outflow impinging on the RV wall. 
¢ Uneven flow distribution through the heat exchanger. 
* Stagnant flow regions beneath the core. 


These flow features can act to oppose the driving buoyancy force from the core, reducing the 
bulk circulation rate. It is difficult to predict the magnitude of these pressure losses and account 
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for them in the system code analysis before carrying out the CFD analysis. Therefore it is not 
unexpected that the CFD model would predict a lower natural circulation flow rate than the system 
code analysis, when making an initial comparison between models. 


The magnitude of the difference between the two models is significant (over 25%) and warrants 
further investigation (along with the unsteadiness in natural circulation flow rate predicted by the 
CFD model). However, taking into account the uncertainty in the prediction of this particular result 
(using either model), and the differences in input data, a sizeable discrepancy is not surprising’. It 
would, therefore, be sensible to align the input data for both models, further refine the CFD model 
and analyse a longer transient period before drawing any firm conclusions. 


If the difference in the natural circulation flow rate between the models is found to be largely due to 
pressure losses caused by the three-dimensional flow features, then these aspects can be updated 
in the system code model to improve its future predictions. 


Reactor Component Temperatures 


CFD (RV Inner Wall) 
e System Code (RV Inner Wall) 


CFD (GV Inner Wall) 

CFD (GV Outer Wall) 

CFD (RV Outer Wall) e System Code (GV Outer Wall) 
° System Code (RV Outer Wall) 


Figure 5.26: Vertical profiles of wall temperature on the outer and inner surfaces of the 
RV and GV. The temperatures are circumferential averages. 


Figure 5.26 compares vertical profiles of wall temperature on the inner and outer surfaces of the 
RV and GV. The system code results that are available are shown as discrete data points, while 
the CFD results are shown as continuous lines. 


Figure 5.26 shows that the inner and outer RV wall temperatures are higher in the CFD results 
than the system code analysis. This observation is consistent with the differences in the predicted 
temperature of the lead pools discussed above. Despite this difference between the CFD and the 
system code results, the temperature profiles in the vertical direction show the same trend, with an 
increase in the temperature of the GV outer wall above the water level (z = 5m). 


It is also noted that both the temporal and spatial discretisation of the models are likely to affect the prediction of the natural 
circulation flow rate. 
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Below the water level, the temperature of the outer surface of the GV is = 100°C in the CFD model 
and ~ 130°C in the system code analysis. The temperature of the outer surface of the GV is not 
specified directly in the model but is generated as part of the boundary condition which is applied 
to the outer surface of the GV to represent the boiling water. Figure 5.26 suggests that there is an 
inconsistency between the boundary conditions applied in the CFD model and the system code 
analysis. An update to this CFD model boundary condition is discussed in Section 6.2. 


Turbulent Prandtl Number 


The choice of RANS turbulence model and especially the turbulent Prandtl number (Pr;) was high- 
lighted in Section 4.4 as one of the main uncertainties in the CFD model setup. Since performing 
CFD modelling with a higher fidelity method such as LES is not practical for this initial study, the 
best way to investigate this issue is via sensitivity studies. 


Referring to Volume 5 (Liquid Metal Thermal Hydraulics), the turbulent Prandtl number is defined 
as the ratio between the turbulent diffusivity of momentum (v;) and the turbulent diffusivity of heat 


(at): 


At 
In their standard implementation, two-equation turbulence models (such as the k-w SST model 
used in this case study) do not separately model the turbulent diffusivity of heat a;. Instead a; is 
calculated from vz: 


Vt 
y= 
Pr; 


Therefore increasing the turbulent Prandtl number has the effect of reducing the turbulent diffusivity 
of heat. The default value for turbulent Prandtl number in ANSYS CFX is 0.85 (in common with 
many CFD codes). As discussed in Volume 5, this value can over estimate the turbulent transport 
of heat in a low molecular Prandtl number fluid where conduction plays a more significant role. 


In the context of this case study, the potential sensitivity to the turbulent Prandtl number is in- 
vestigated by performing three ‘steady-state’ simulations with turbulent Prandtl numbers of 0.85, 
2.0 and 4.0 respectively. All other model inputs (including the RANS turbulence model and wall 
functions) were identical. 


From the equations above, it might be assumed that changing Pr; would not have a significant 
effect on the predicted momentum field. However, in the case of natural circulation, where the 
momentum of the flow is driven by temperature differences, changing the turbulent Prandtl number 
is likely to affect both the flow and temperature distributions within the model. 


Figure 5.27 shows contours of temperature on a vertical slice through the lead pools, for turbulent 
Prandtl numbers of 0.85, 2.0 and 4.0. These contours were extracted from ‘steady-state’ solutions 
and do not fully resolve the unsteady flow features. As shown in Section 5.3, the unsteady flow 
features have a significant effect on the velocity and temperature fields in some areas of the lead 
pools and it is therefore difficult to attribute differences in the temperature contours to the turbulent 
Prandtl number alone. 
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(a) Pr, = 0.85 (b) Pr, = 2.0 (c) Pr; = 4.0 
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Figure 5.27: Contours of temperature on a vertical slice through the lead pools, for tur- 
bulent Prandtl numbers of (a) 0.85, (b) 2.0 and (c) 4.0. 


Nevertheless, there are noticeable differences in the temperature contours in the pools. In the hot 
pool for example (where the transient variations in the pool were found to be low), the hot plume 
near the reactor centreline mixes and cools less in the analysis with a turbulent Prandtl number of 
4.0 than at a turbulent Prandtl number of 0.85. It can also be seen that the temperature variation 
as a whole in the hot pool is reduced in the analysis with Pr; = 0.85 compared with the other 
two cases, suggesting more dissipation in the thermal field (as would be expected). In the lower 
cold pool, there also appears to be differences in the stratified flow region underneath the core 
and its interaction with the density-driven current. However, given the transient variation in the flow 
patterns in this region, it is difficult to draw any firm conclusions. 


For a closer examination of the hot pool, Figure 5.28 shows contours of temperature and velocity 
on a horizontal slice at a height of z =4m, for turbulent Prandtl numbers of 0.85, 2.0 and 4.0. Near 
the centreline, the temperature of the plume is predicted to be = 2°C higher at a turbulent Prandtl 
number of 4.0 than at a turbulent Prandtl number of 0.85. The temperature gradient between the 
centreline of the hot pool and the bulk of the pool is also predicted to be steeper when Pr; = 4.0 
than when Pr; = 0.85. This observation is consistent with the decrease in thermal diffusivity that 
occurs when the turbulent Prandtl number is increased. 
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Figure 5.28: Contours of temperature (top) and velocity (bottom) on a horizontal slice through the hot pool at z = 4m, 
for turbulent Prandtl numbers of (a) 0.85, (b) 2.0 and (c) 4.0. 


In the bulk of the pools, the velocities are relatively low (~ 0.1 m/s) and the turbulence levels are 
small. As a result, the turbulent transport of heat (and the turbulent Prandtl number) has a relatively 
small effect on the majority of the temperature field. The contours in Figures 5.27 and 5.28 confirm 
that, even in the hot pool, uncertainty in the turbulent Prandtl number has a relatively small effect 
on the temperature of the flow through the heat exchanger into the upper cold pool. 


These contour plots are supported by the bulk (volume-average) temperature of the lead pools and 
the RV surface temperatures, which are found to be less than +3°C different from the baseline case 
across the range of Pr; analysed. It follows that the uncertainty in the turbulent Prandlt number is 
likely to have a small effect on the overall predictions of primary circuit temperature (especially 
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when compared with the uncertainty in the model boundary conditions). However, the possibility of 
it having a larger local effect in specific areas can not be ruled out at this stage. 


Pr; Circulation Rate (kg/s) Difference (%) 
0.85 1854 

2.0 1729 6.7 

4.0 1625 12.4 


Table 5.2: Circulation rate of primary coolant for turbulent Prandtl numbers of 0.85, 2.0 and 4.0. 


Examining other aspects of the results, it was noticed that the predicted natural circulation flow 
rate (shown in Table 5.2) appears sensitive to the value of turbulent Prandtl number, with a higher 
value being predicted for lower values of Pr;. However, considering the lack of convergence in the 
prediction of this parameter in the steady-state analysis, it would be of benefit to consider a study 
using transient analyses to investigate this further. 


It should also be noted that this sensitivity study does not, in itself, shed any light on what the best 
value of turbulent Prandtl number is to use for this model. As discussed in Volume 5, Section 3.3.2, 
the manipulation of Pr; in this way is a compromise and possibly not one that works especially 
well for complex buoyancy driven flow. If further study reveals that it does have a significant effect 
on the results of importance then, in the absence of more sophisticated methods, this may be an 
uncertainty that persists in future CFD models. 


Pool Free Surface Boundary Condition 


As the lead pools have a large surface area, the heat loss from the surface of the pools may make 
a significant contribution to the overall heat balance for the PHRS. Heat will be lost from the surface 
by both convection to the cover gas and radiation to the structures above the surface of the lead 
pools. Details of the cover gas and above core structure were not available to provide inputs for the 
analysis, so a simple thermal boundary condition was adopted for the top surface of the lead in the 
CFD model. However, this boundary condition is subject to considerable uncertainty. 


Interrogation of the baseline model predictions reveals that heat transfer from the surface of the 
pools comprises over 10% of the total thermal energy lost from the model. Even at lower lead 
pool temperatures (more representative of normal operation), a considerable parasitic heat loss 
through this surface would still occur if the current boundary conditions are maintained. Although 
the design of this area of the reactor is not complete, it is unlikely that this magnitude of heat loss is 
realistic. Furthermore, the flow patterns in the hot and upper cold pools suggest that this boundary 
conditions is having an affect of the phenomenology of the flow, potentially reducing stratification 
in a way that may compromise the objectives of this study. 


To help quantify the effect of the uncertainty in the pool surface heat transfer boundary condition, 
the CFD model was rerun with the thermal boundary condition on the lead surface replaced with 
an adiabatic boundary condition. This represents the bounding case where no heat is lost from the 
free surface and all the heat is removed through the outer surface of the GV. 


Figures 5.29 and 5.30 show comparisons of temperature contours and velocity vectors on a vertical 
plane through the lead pools, with the original and adiabatic boundary conditions applied to the free 
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Figure 5.29: Comparison of primary circuit temperature distribution with (a) the original 
and (b) an adiabatic free surface boundary condition. Note the different scales 
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Figure 5.30: Comparison of primary circuit flow patterns with (a) the original and (b) an 
adiabatic free surface boundary condition. 
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surface of the lead. As would be expected, the removal of the heat loss path through the lead free 
surface leads to a global increase in the bulk lead temperature (of around 40 °C). 


The temperature distribution and, in some areas, the flow patterns in the hot pool, are also locally 
affected by the adiabatic free surface boundary condition. In both cases a buoyancy-driven plume 
rises from the centre of the core to the top of the pool. With the original boundary condition, this 
hot flow is cooled as it travels radially outward along the free surface and starts to mix. In the case 
with the adiabatic boundary condition, this warmer flow persists for longer, until it gets close to the 
cooler core barrel surface. The variation in temperatures with height in the hot pool away from the 
central plume is higher in the adiabatic case than in the original model and the hot pool exhibits 
reduced mixing, lower flow velocities and increased stratification when the surface is adiabatic. In 
the upper cold pool the differences are more subtle, but signs of reduced mixing near the top of the 
pool can still be observed in the predictions with the adiabatic surface boundary condition. 


Despite these local differences, there is still no evidence of major additional areas of stratification 
that may be problematic to PHRS performance. The flow patterns and temperature distribution 
in the lower cold pool appear largely unaffected (within the range that would be expected from 
a solution with unsteady flow features) and notable aspects of the flow pattern, such as the cool 
current of flow down the RV wall, are still present. The magnitude of the temperature differences 
between the pools and range of temperatures within the pools is very similar between the two 
analyses. 


It is concluded that the qualitative effect on the results of changing between the two free surface 
heat transfer boundary conditions is reasonably small. However, the magnitude of the affect on the 
(quantitative) bulk primary circuit lead temperatures is more significant, and determining a more 
realistic boundary condition should be given a high priority. 
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Application of the Results and Next Steps 


The aim of this case study has been to demonstrate how to use thermal hydraulic modelling (in 
this case CFD) to contribute to early stage reactor design activities in a representative ‘real world’ 
situation. To achieve this objective, the emphasis has been placed on demonstrating the planning, 
execution and interpretation of an ambitious initial model, rather than providing a demonstration of 
how a refined outcome can be achieved with a more mature model. 


Ideally a CFD model of this sort would be continuously developed and refined alongside the reactor 
design. However, in order to justify the continuous development of such a model to stakeholders, 
the model needs to provide continuous value to the reactor design project at each stage. To max- 
imise this value, the modelling engineer should have a clear understanding of the strengths and 
weaknesses of the model and a vision for moving it forward. 


To illustrate this process, the application of the current model results and the next stages in the 
model development are described in the following sections: 


* Section 6.1 describes the conclusions from the current model, both in terms of input to the 
Ongoing reactor design process and what has been learned regarding the use of CFD mod- 
elling in this case. 


* Section 6.2 describes improvements that could be made to the model itself to improve its 
accuracy. 


* Section 6.3 discusses the next steps that could be taken in the overall model development 
process. This section includes simple extensions to the model itself to expand its capability 
and also additional modelling that could be carried out to provide further insight. 


Conclusions from the Analysis 


Conclusions for Reactor Design and Safety 


The CFD modelling undertaken has provided considerable insight into the three-dimensional flow 
of lead inside the primary circuit of a LFR under decay heat removal conditions. The flow visual- 
isation provided by a CFD model, used in combination with a knowledge of fluid mechanics and 
heat transfer, can significantly improve understanding, even if the model simply confirms previous 
findings. 


With a plethora of potentially available results, it is easy (and tempting) to extract and utilise a wide 
range of results from the CFD model. Therefore, it is important that the developer of the CFD model 
provides guidance to engineers of other disciplines regarding what can and cannot be reasonably 
concluded from an early stage model of this type. From the examination of the CFD model results 
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in Section 5, taking into account the limitations of the current model, it is sensible to primarily draw 


qualitative conclusions rather than quantitative conclusions as follows: 


The analysis confirms the operation of the PHRS under passive heat removal conditions 
is in line with the intention of the design. The overall variation in lead temperatures and 
flow velocities are largely in line with the results of previous calculations (with a few notable 
difference discussed further in Section 6.1.2). 


When the PHRS heat sink is full of water, especially boiling water, then the rate of heat 
transfer from the primary circuit is dominated by the RV-GV gap and conduction through the 
steel vessels. That is to say that the HTCs between the lead and the RV and between the 
GV and the water are much higher than the effective HTC across the gap. Hence, the vast 
majority of the temperature difference between the primary circuit lead and the environment 
is across the RV, the gap and the GV and changes to their heat transfer characteristics will 
have the largest effect on reactor coolant and vessel temperatures. 


The majority (~ 90%) of the heat transfer to the PHRS heat sink, across the RV-GV gap 
under passive heat removal conditions occurs by thermal radiation, with natural convection 
of the gas within the gap playing a much smaller role. Having a good (or bounded) knowl- 
edge of the emissivity of these surfaces (in the reactor and in any future test facilities) will 
be important to support PHRS safety arguments. It is noted that thermal radiation will also 
be important to the heat transfer on the outside of the GV above the water level and so 
knowledge of the emissivity of the external surface of the GV is also likely to be important. 


The current design will almost certainly experience thermal stratification of the lead below 
the core when the PHRS heat sink contains boiling water. Whilst there is no indication from 
the current CFD analysis that this is particularly detrimental to the successful operation of the 
PHRS, it is a potential contributor to the lower natural circulation rate observed in the CFD 
analysis compared with the system code analysis. Structural components passing through 
this stratified region will see corresponding temperature gradients. 


Other areas of possible stratification within the primary circuit have been identified, but com- 
prise comparatively small temperature variations (~ 10 °C) from the bulk temperature of their 
respective pools. 


Below the water level, a density-driven current of cooler lead is predicted to flow down the RV 
wall, whilst the flow patterns in the middle of the cold pools, away from walls, are dominated 
by gently swirling flow. The presence of a region of cooler lead (approximately 20°C cooler 
than the lead in the majority of the pool) will affect the heat transfer between the primary 
circuit fluid and the RV wall and may need to be taken into account in the system code model 
in the future. 


The results of the CFD modelling suggest that the bulk natural circulation flow rate in the 
primary circuit lead may be both unsteady and difficult to quantitatively predict with a high 
level of confidence. However, the model predictions support that the PHRS design is capable 
of generating a naturally circulating flow in the lead coolant such that the bulk temperature 
variation throughout the primary circuit is only a few 10s of degrees (outside of the stratified 
region at the bottom of the reactor). 


Examination of the flow through the heat exchanger and the core identified that the predic- 
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tions are likely influenced by model simplifications (explored further in Section 6.2). Despite 
these issues, it is likely that the flow through these components will not be evenly distributed 
under naturally circulation conditions. There is no indication that this will be detrimental to 
the performance of the PHRS, but will make the pressure losses in the primary circuit (and 
therefore the natural circulation flow rate) harder to estimate in simpler calculations. 


* The comparison between the CFD and system code modelling highlighted the sensitivity of 
the bulk lead temperature predictions of both models (and likely the real system performance) 
to: 

— The surface area of the RV and GV. 

— Assumptions regarding the rate of boil down of the PHRS heat sink. 
— The variation in core decay heat during the fault. 

— Heat losses from the upper pool surfaces. 


These findings may be relevant to the consideration of other faults and/or the evolution of the 
design of the reactor and PHRS. 


Conclusions for Modelling 


In addition to providing insight into reactor operation under PHRS conditions, this model was also 
intended to explore the potential benefits of using CFD in conjunction with system code modelling 
in the future. It is, therefore, also of value to consider what conclusions can be drawn from a 
modelling perspective to help inform future plans: 


* Overall it is concluded that CFD is a useful addition to the predictive tool set to support reactor 
design and safety substantiation. The model has been successful in fulfilling its objectives 
with regard to supplementing the system code predictions in specifically identified areas and 
provided further evidence of the success of the PHRS design. 


* The comparison between the steady-state and transient calculations concluded that the over- 
all flow patterns were reasonably well captured by a steady-state approach. However, the 
poor convergence of the steady-state momentum field suggests that there are unsteady ther- 
mal hydraulic phenomena present and it would be difficult to use a steady-state CFD model 
solution with confidence. It is therefore important that future CFD studies include provision 
for performing transient analyses. 


« The timescale and computational resources available to the engineer were not sufficient in 
this case to gain maximum value from this CFD model. Only two iterations of the mesh were 
possible in the timescales, the duration of the transient analysis was shorter than would be 
preferable and a number of model improvements and sensitivity studies remain outstanding 
(see Section 6.2). However, this is not an unusual situation for a ‘first’ model of a specific 
application. A break point at this stage (whilst frustrating for the modelling engineer) allows 
other technical areas of the project to both provide input into and start to take value from the 
work and allows engineering managers to consider the potential benefits of the modelling 
and incorporate it into future plans. 


¢ The sensitivity studies found that varying the turbulent Prandtl number made very little dif- 
ference in the overall system temperatures most important to PHRS performance. However, 
more variation was observed in the predicted natural circulation flow rate and in the detailed 
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temperature distribution within the pools. Therefore it cannot be concluded that the correct 
prediction of the turbulent dissipation of heat is unimportant in this case and it may be more 
important under conditions such as normal operation, where the flow will be more turbulent. 
It is therefore recommended that future CFD modelling continues to consider this an area 
of uncertainty and appropriate studies are performed. Additionally, there is likely to be ben- 
efit in continuing to maintain awareness of ongoing research that is investigating RANS and 
URANS modelling of heat transfer in liquid metals to facilitate continuous improvement of this 
aspect of the model. 


¢ The sensitivity studies show that the heat transfer from the free surface of the reactor pools 
can have a significant effect on the overall bulk lead temperatures and a small, but noticeable, 
effect on the mixing in the hot and upper cold pools. Further analysis to investigate the 
presence and/or influence of thermal stratification should, therefore, attempt to include a 
more realistic boundary condition on this surface. 


Improvements to the Current Model 


Section 5 highlighted regions of the CFD model where the predictions were considered to have 
been sub-optimal due to the boundary conditions, mesh resolution or geometrical simplification. A 
sensible next step in progressing the CFD model is to address as many of these issues as possible. 
However, it is important to be mindful of the investment in time required to achieve improvements 
and weigh this against the potential value gained. For example, in this case there will inevitably 
be considerable uncertainty in the quantitative predictions of the model due to missing input and 
validation data. Therefore, it is sensible to hold back from carrying out too much fine tuning of the 
model at this early stage. Relatively simple improvements and activities that will help to focus future 
studies should be the priority. 


Improvements to the Mesh 


A number of areas in the mesh were identified for improvement in Section 5, all of the which 
require increased resolution. Hence, the required improvements will all lead to an increase in the 
total number of cells in the mesh and therefore an increase in the required computing time for the 
model. It should be cautioned that if the mesh becomes too large and cannot be analysed in a 
reasonable time with the available computing resources, then the value of the model to the project 
will be compromised. It is also reiterated that the pursuit of a proven mesh independent solution for 
this model would be challenging and may not be justified at an early stage of the reactor design. In 
this case, the following ‘short-term’ mesh improvements are recommended: 


¢ Refinement of the boundary layer mesh on the RV wall in the lower cold pool (and potentially 
the upper cold pool as well) is likely to result in increased accuracy of the heat transfer 
calculations. 


¢ Improvements to the cell size and distribution near the bottom of the upper cold pool would 
provide greater confidence in the prediction of stratification in this region and flow distribution 
through the heat exchanger. 


* Refinement of the mesh in the region below the core (in combination with improvements to 
the representation of the geometry in this region) would improve confidence in the prediction 
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of the stratification and interactions with the density-driven current in this area. 


* Increasing the number of cells through the thickness of the RV is recommended to avoid 
a very large discrepancy in the mesh size at the fluid-solid interface. Increases in mesh 
resolution in a vertical direction should also be considered in areas where there are higher 
vertical temperature gradients in the adjacent lead. 


¢ Although not identified as a particular area of concern, it is judged that the mesh in the 
boundary layer adjacent to the sides of the core barrel is also likely to need some refinement 
(although not as much as that on the RV wall). 


Improvements to the Representation of the Geometry 


Although in most cases the simplifications to the geometry in this model have been necessary due 
to a lack of available information, there are some areas where it has been identified that these 
simplifications may have a significant effect on the results. 


The most notable example of this is the geometry of structures at the entrance to or exit from 
the core. At the current time, this information is not available but it will impact the predictions of the 
overall natural circulation flow rate and the velocity and temperature distribution below, through and 
above the core. There may be little that can be done about this in the short-term but it is highlighted 
here as an area of high priority for the future. 


The heat exchangers are currently represented by a single block of porous material. Although 
the detailed geometry of this component is not available, it is understood that it will comprise 
a number of discrete channels through which the lead will flow separated by thick plates which 
contain internal passages for the secondary coolant. Interrogation of the CFD results revealed that, 
despite the increased resistance in the vertical and circumferential directions, the flow through the 
heat exchanger redistributes in a manner that would not be plausible in reality. Inserting a number 
of thin walls within the porous media would limit the ability of the flow to move laterally without 
requiring significant changes to the mesh. 


Improvements to the Model Boundary Conditions 


The upper surface of the lead pools is the only location where thermal energy is removed from the 
model, other than the PHRS heat sink. The current representation of this heat loss is simple and 
is likely to result in an over-estimate. The sensitivity study investigating this boundary condition 
shows that it has some influence on the predictions of temperature distribution and the level of 
convection driven mixing in both the hot pool and the upper cold pool and a larger influence on 
the overall predicted lead temperature. It is therefore worth subjecting this boundary condition to 
further scrutiny, improving the estimate of the temperature to which this surface is radiating and 
potentially including a provision for convective heat transfer. 


The boundary conditions applied to the outer surface of the GV to represent the PHRS heat sink 
are described in Section 4.1.7. Although these have been subject to some consideration, a number 
of simple improvements could be implemented to improve the model. The precise value of the HTC 
below the water level is unlikely to be important, as it is much larger than the effective HTC across 
the RV-GV gap. However, it could be improved by incorporating a temperature dependent HTC and 
the variation in the water boiling point with hydrostatic pressure. 
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The limitations that arise from using a sector model of the reactor (with periodic boundary condi- 
tions applied to the sides of the sector) has been discussed in context throughout Section 5. How- 
ever, there is little that can be done to improve this modelling approach without greatly increasing 
the total number of cells in the mesh and the computational resources needed to solve the model. 
Whilst it is unlikely that this could be addressed in the short-term, longer term suggestions are 
made in Section 6.3. 


Additional Analyses and Future Modelling 


The initial analysis has highlighted the utility of using CFD analysis to support calculations from 
system code analysis. Therefore, the CFD model is likely to now form part of the wider analysis 
toolset used for continuing to develop and substantiate the reactor design. Having invested the time 
into creating and refining a new CFD model of the reactor it is likely that the project team will wish 
to put it to further use. Describing all of the possible applications of CFD modelling to LFR reactor 
development is beyond the scope of this document. However, a number of ideas are presented 
below as examples of what might be done next. 


Simple Extensions Using the Current Model 


The model described in this document was developed with the broad objective of investigating the 
influence of three-dimensional flow and stratification on PHRS operation under a specific set of fault 
conditions. The analyses completed to date have met this objective, but only at a single point in 
time after the initiation of the fault. Previous publications (See Section 3.5) suggest that in the early 
stages of the fault, the CFD model is less likely to be successful in producing accurate predictions 
(and a ‘snapshot’ modelling approach would not be appropriate). However, at later points in time, 
when the water level has dropped a little further, it would be reasonably straightforward to modify 
the existing model so that it can be used to investigate the changes in the flow and temperature 
distributions in the primary circuit. 


As the water level drops and more of the GV is exposed, the boundary condition applied to the 
external GV wall above the water level will have a more significant effect on the temperature dis- 
tributions in the PHRS. Under these conditions, the current formulation of the external boundary 
condition on the outer surface of the GV (an applied HTC) is likely to be inaccurate. As the water 
level drops it will therefore become necessary to incorporate a more detailed representation of the 
PHRS heat sink cavity, either explicitly in the model or as a more complex set of boundary condi- 
tions. It follows that it is likely to be more complicated to investigate the points in time when most 
(or all) of the water has boiled away, than when the water is, for example, 1m lower than the level 
already analysed in this case study. 


The only sensitivity study carried out with the CFD model in this case study relating to model set-up 
has been a steady-state investigation of turbulent Prandtl number. The current model could also 
be easily used to investigate the sensitivity of the predictions to a range of input parameters, that 
would provide valuable information for both this and a range of future CFD models relating to this 
reactor. 
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Examples of such studies include: 


¢ Mesh resolution in key areas. 

¢ Time-step resolution for the unsteady analysis. 

* Pressure loss coefficients of various reactor components. 
¢ Turbulence modelling approach. 

¢ Wall modelling approach. 


Looking at the list above (which only includes a small number of examples) it is easy to imagine 
that the total number analyses needed to investigate model sensitivity could become large. There- 
fore, it is important to consider the value of what might be learned in each case. For example 
(unless needed for design purposes), it might be better to wait until further information about im- 
portant reactor components becomes available before investigating model sensitivity to pressure 
loss coefficients, so that the study can be better focused. It also may be more efficient to carry out 
preliminary studies into component losses using the system code model’. 


A different issue is well illustrated by the turbulent Prandtl number sensitivity described in Sec- 
tion 5.4.2. Although it is useful to be aware that the model is sensitive to this parameter, the study 
does not enable any conclusions to be drawn regarding the ‘best’ value to use. As the use of a 
modified turbulent Prandtl number in a two-equation turbulence model is a compromise in itself, 
even a more extensive study of this kind (investigating a wider range of values for Pr) would not 
be able to conclusively identify the best value to use (and would not necessarily bound reality). A 
more extensive study is, therefore, not likely to be justified in the absence of any validation data. 


In contrast, investigations relating to turbulence and wall modelling could be carried out using a 
much simpler sub-model and compared with LES predictions in a similar manner to that demon- 
strated by Study A (which provides a good starting point). The findings of this type of study would 
be useful to a variety of future CFD models and there is no technical reason not to progress these 
at the earliest opportunity. Similarly studies examining the sensitivity to mesh and time-step res- 
olution are both a high priority and have the potential to inform a variety of future CFD modelling 
tasks. 


Volume 4 (Confidence and Uncertainty) provides further information on the use of sensitivity stud- 
ies for improving confidence in model predictions. 


Expanding the Model Capability 


In addition to describing the important results that the initial CFD analysis aimed to predict, Sec- 
tion 3.1 also summarised a number of other results which were considered likely future targets for 
this model. In all cases, achieving these results requires some modifications to be made to the 
model. However, the difficulty of making these modifications varies considerably. 


For example, predictions of the three-dimensional temperature and velocity distribution during nor- 
mal, full-power reactor operation will be significantly different to those predicted under natural circu- 
lation conditions. In addition to changes in the model input data (such as core heat load), analysis 
of normal operation would require changes to made to the representation of both the pump and 
the heat exchanger, both of which would be operational. 


1 CFD sub-models could be used to better define component pressure loss once sufficiently detailed geometry is available. 
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In the absence of the system code model, the changes to the heat exchanger required to predict the 
distributions of heat flux and temperature would be quite challenging. However, if the intention was 
again to use CFD to supplement system code analysis results (rather than provide an independent 
check), then heat flux information could be directly taken from system code model results and 
applied to the porous media in the CFD model, significantly simplifying the process. 


Incorporating a detailed model of pump impeller and guide vanes into the current CFD model would 
be a significant undertaking and greatly increase its complexity. However, if the pump performance 
is well defined, a simple approach to modelling the operational pump would be to apply a momen- 
tum source to the current block of porous media. This approach would be unable to capture some 
of the flow details, but could include swirl generated by the pump. 


Additional CFD Modelling 


Throughout this document various additional models/sub-models have been suggested to inves- 
tigate specific aspects of the reactor thermal hydraulic performance under PHRS operation. It is 
not essential that all of these additional models are created, but it is useful to collate the ideas 
generated as part of the work to date, so that they can be called upon if the need arises. Three 
examples of the types of models that may be considered are included below. 


Sub-models to Investigate Convective Heat Transfer Modelling in Lead: The use of a two- 
equation RANS turbulence model and wall functions for the modelling of boundary layers contain 
approximations when applied to any fluid. As previously discussed, when being used to model a 
low Prandtl number fluid their applicability carries increased uncertainty. In this case, many of the 
surfaces in contact with the primary circuit where heat transfer is important are vertical (or curved) 
and will have a buoyancy influenced boundary layer, with the ‘free stream’ flow being extremely 
slow. Study A provides an example of a simpler geometry that includes flow over a heated vertical 
section and provides a good starting point. However, more application specific, high-fidelity (e.g. 
LES) sub-models examining, for example, the flow down sections of the RV wall would be extremely 
helpful in determining the applicability of the methods employed in this case study. 


Sub-models to Investigate Phenomena in the Lead Pools: Developing stand-alone models 
of each of the reactor pools to investigate natural circulation conditions would have been very 
challenging at the outset of this work due to the two-way coupling between the flow and heat 
transfer around the whole system. However, the ‘whole primary circuit’ model, developed within 
this case study, is now available to provide representative boundary conditions to smaller models. 
Sub-models of each reactor pool have the following advantages: 


* Each sub-model can be created with a more refined mesh (whilst still being reasonably quick 
to analyse) to incorporate more detailed geometry (when it becomes available) and to per- 
form mesh sensitivity studies. 


¢ For regions such as the hot pool, where it is considered more likely that the periodic boundary 
conditions may be influencing the model predictions, a full 360° model could be built. 


* Other types of sensitivity study, to investigate user assumptions and modelling choices, could 
be performed more efficiently using a smaller model. 
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Modelling to Investigate the PHRS Heat Sink: The modelling performed to date has used 
boundary conditions to represent the PHRS heat sink. When the majority of the GV surface is 
covered with boiling water, this approach is likely to be reasonable (i.e. the HTC is so high that the 
model predictions are not very sensitive to its exact value). However, as the water boils down and 
the GV is exposed, the convective and radiation heat transfer is much hard to represent in this way. 


A separate model of the PHRS heat sink including the full cavity geometry, solid structures and air 
chimneys, with representative boundary conditions on the RV inner surface determined from the 
reactor primary circuit model (developed by this case study), would give considerable insight into 
the cooling above the water surface and, later in fault, when all the water has boiled away and the 
GV is cooled by air. Additionally, it can be used to provide insight into the temperature distribution 
and pressure losses within the PHRS heat sink system itself to give a more complete picture of the 
system performance. 
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8 Nomenclature 


Latin Symbols 


— 


Area, m? 

Atwood number (At = (¢1 — p2)/(p1 + p2)) 

Biot number (Bi = hL/k,) 

Specific heat at constant pressure or volume, Jkg~! K~+ 

Diameter (Dp = 4Acs/Pes for hydraulic diameter), m 

Darcy-Weisbach friction factor 

Fourier number (Fo = a,t/L?) 

Grashof number (Gr = gL?Ap/v?p = gL°BAT/v, using the Boussinesq approxi- 
mation Ap/p + —BAT, where AT is often taken as Ty — T5,00) 
Acceleration due to gravity, ms~? 

Specific enthalpy, J kg~+, Heat Transfer Coefficient (HTC), Wm? K~? or height, m 
Radiative intensity, W m~? sr~! or Wm? sr~! wm? for a spectral density, where sr 
(steradian) is solid angle 

Radiosity, Wm~? 

Thermal conductivity, W m7! K~! 

Length or wall thickness, m 

Molar mass of a species, kg kmol~* 

Mach number (Ma = U/a, where a is the speed of sound) 

Refractive index 

Nusselt Number (Nu = hL/kg) 

Perimeter, m 

Pressure (P; = static pressure, Pr = total pressure), Nm~? or Pa 

Péclet number (Pe = RePr) 

Prandtl number (Pr = cpu/k¢) 

Heat flux (rate of heat transfer per unit area, gq = Q/A), Wm~? 

Rate of heat transfer, W 

Radius, m 

Gas constant (for a particular gas, R = R/M), Jkg~! K~! 

Universal gas constant, 8314.5 Jkmol~! K+ 

Thermal resistance, K W~4 

Rayleigh number (Ra = Gr Pr) 

Reynolds number (Re = pUL/, or for an internal flow Re = WD»,/A-sL) 
Richardson number (Ri = Gr/Re’) 

Strouhal number (Sr = fL/U, where f is frequency) 
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Stanton number (St = Nu/RePr) 

Time, s 

Temperature (7, = static temperature, 77 = total temperature), K 
Wall friction velocity (u; = \/Tw/p),ms~+ 

Velocity, ms~! or thermal transmittance, W m~? K~4 
Specific volume, m? kg~! 

Volume, m? 

Mass flow rate, kgs~! 

Wall distance, m 


Non-dimensional wall distance (y* = yu, /v) 


Greek Symbols 


yn 2 DBD RQ 


Qn 8S FY BR 


4 


2-1 


Thermal diffusivity (a = k/pcp), m*s~ 
Volumetric thermal expansion coefficient (8 = —(1/p)(0e/OT)), K~? 
Ratio of specific heats (y = cp/cv) 

Emissivity or surface roughness height, m 

Absorption coefficient, m+ 

Wavelength, m 


Viscosity, kg m~!s~1 


Kinematic viscosity and momentum diffusivity (v = u/p), m?s~? 


Density, kg m~? 
Stefan Boltzmann constant, 5.67 x 10-? Wm-2 K~4 
Shear stress, Nm~? 


Porosity or void fraction 


Subscripts and Modifications 


b 
cs 
f 
i 
+ 
t 


Bulk (mass-averaged) quantity 
Cross-sectional quantity 

Quantity relating to a fluid 

Quantity relating to a particular species 
Total (stagnation) quantity 

Turbulent quantity 

Static quantity or quantity relating to a solid 
Quantity relating to a wall or surface 
Quantity far from a wall or in free-stream 
Average quantity 

Molar quantity 

Varying quantity 
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9 Abbreviations 


AMR 
CAD 
CFD 
CFL 
CHT 
DNS 
FCH 
GV 
HTC 
LES 
LFR 
LMFR 
LWR 
NPP 
PHRS 
PIRT 
RANS 
RMS 
RV 
S2S 
SFR 
SST 
URANS 


Advanced Modular Reactor 
Computer Aided Design 
Computational Fluid Dynamics 
Courant-Friedrichs Lewy 

Conjugate Heat Transfer 

Direct Numerical Simulation 

First Cell Height 

Guard Vessel 

Heat Transfer Coefficient 

Large Eddy Simulation 

Lead-cooled Fast Reactor 

Liquid Metal-cooled Fast Reactor 
Light Water Reactor 

Nuclear Power Plant 

Passive Decay Heat Removal System 
Phenomena Identification and Ranking Table 
Reynolds-Averaged Navier-Stokes 
Root Mean Square 

Reactor Vessel 

Surface-to-Surface 

Sodium-cooled Fast Reactor 

Shear Stress Transport 

Unsteady Reynolds-Averaged Navier-Stokes 
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